{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "This notebook looks at printing a vertical profile in EFDC.\n",
    "\n",
    "\n",
    " \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "import pandas as pd\n",
    "import netCDF4 as nc\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/fearghalodonncha/opt/anaconda3/lib/python3.8/site-packages/ipykernel/ipkernel.py:287: DeprecationWarning: `should_run_async` will not call `transform_cell` automatically in the future. Please pass the result to `transformed_cell` argument and any exception that happen during thetransform in `preprocessing_exc_tuple` in IPython 7.17 and above.\n",
      "  and should_run_async(code)\n"
     ]
    }
   ],
   "source": [
    "# generally in EFDC our domain is bounded by land cells at the outer edges; drop these from viz\n",
    "land_bounded_cells = 2  \n",
    "def getcoordinates(ncdataset):\n",
    "    x = ncdataset.variables['X'][:]\n",
    "    y = ncdataset.variables['Y'][:]\n",
    "    x = x[land_bounded_cells:-land_bounded_cells]\n",
    "    y = y[land_bounded_cells:-land_bounded_cells]\n",
    "    time = ncdataset.variables['Time'][:]\n",
    "    return x, y, time\n",
    "\n",
    "def getvariablevalues(ncdataset, variable):\n",
    "    var0 = ncdataset.variables[variable][:]\n",
    "    if (len(var0.shape)==1):\n",
    "        return var0.data\n",
    "    var1 = var0[0]\n",
    "    if (len(var1.shape) == 2):\n",
    "        var1 = var1[land_bounded_cells:-land_bounded_cells, land_bounded_cells:-land_bounded_cells]\n",
    "    if (len(var1.shape) == 3):\n",
    "        var1 = var1[:, land_bounded_cells:-land_bounded_cells, land_bounded_cells:-land_bounded_cells]   \n",
    "    return var1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-3-268344256d76>:4: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  x = ncdataset.variables['X'][:]\n"
     ]
    }
   ],
   "source": [
    "data_directory = \"../../SampleModels/CobscookBayModel/\"\n",
    "fname = data_directory + \"efdcout_20190102000000.nc\"\n",
    "ncdataset = nc.Dataset(fname)\n",
    "x, y, time = getcoordinates(ncdataset)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Current oversight is that water depth are not included in NetCDF files\n",
    "Read these instead from DXDY.INP <br>\n",
    "Clunky way of doing it but works for now\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Current oversight is that water depth are not included in the netcdf file\n",
    "## We need to read from DXDY.INP file\n",
    "dxdy = pd.read_csv(data_directory + \"DXDY.INP\", skiprows=4,\n",
    "                   names = [\"I\", \"J\", \"dx\", \"dy\", \"dpth\", \"elev\", \"zruff\", \"veg\"], \n",
    "                   delim_whitespace = True)\n",
    "## Create the 2D array based on I and J coordinates\n",
    "## remember to account for fact we:\n",
    "## 1) stripped two cells along the edges (-2)\n",
    "## 2) Python is zero-based (-1)\n",
    "depth = np.zeros([len(y), len(x)])\n",
    "for il in range(0, len(dxdy)):\n",
    "    c = dxdy[\"I\"].iloc[il] - 2 - 1\n",
    "    r = dxdy[\"J\"].iloc[il] - 2 - 1\n",
    "    depth[r, c] = dxdy['dpth'].iloc[il]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Total depth and plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-3-268344256d76>:12: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  var0 = ncdataset.variables[variable][:]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<matplotlib.colorbar.Colorbar at 0x7fb66e7b98b0>"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXEAAAEKCAYAAADkYmWmAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAA8cUlEQVR4nO29f5wddXX//zwEdiMkJKRZScoSAiEJRT4UIYKVUhOQllgEy9cf8NGaIjZfbURQqZDSCn588AHUYvEL1W8ECiiCK2JFaNQY2aIf5UeCEAgEQiSkK4EQ+bnYbAw5nz9mZjM7O7/vzJ2Zvef5eNzH3pk7P86du/c15573eZ8jqophGIbRTHar2gDDMAwjPybihmEYDcZE3DAMo8GYiBuGYTQYE3HDMIwGYyJuGIbRYEzEDcMwKkJExonIr0TkDnf5CBG5R0QeFJFVInJ00jFMxA3DMKrjHOAx3/IXgM+p6hHAZ93lWEzEDcMwKkBEeoG/BK7xrVZgb/f5JOCZpOPsXrxp7WPy5Ml68MEHV21GLK+99hp77bVX1WbEUoWN6x/aNGJ59h/PiN3ermPr1N0+KNfG1atXb1XVnlaOcdz88friCztTbbv24d+vBbb5Vi1T1WW+5X8BPgNM9K07F/iRiHwJx8l+W9J5Gi3i++67L6tWrarajFj6+/uZP39+1WbEUpWNC6cvGX6+fNXVsdvadWydutsH5dooIk+3eowXX9jJd++cmmrbQ2Zs3qaq8yJsORnYoqqrRWS+76WPAZ9U1e+KyPuAa4F3xJ3HwilGZSzfHC/chjGGORY4RUQ2ArcAx4vIN4FFwG3uNt8BbGDTqDeekC+cvmSEZ24YYxlVXaqqvao6Ezgd+KmqfhAnBv52d7PjgfVJx2p0OMUw2kHw5mK/IIwS+VvgShHZHSeevjhpBxNxo3KWb766tl74iLi9ibdRAqraD/S7z38OHJVlfwunGEYEdb2xGIYfE3Gjcpoilha3N+qIibhRK+okkhY+MZqAxcSNtuKJtD8rJWobP+ff/N5yDQuhTjcUw4jCRNyoBBNIwygGC6cYpRAVP64iRJE3lh1lq92AjDphnrjRdtKkFAYFtL+/P9e5WhFcE2ujCZiIG5UQFxOvekAx6w3GMKrEwilGqeQRxKI8YPOkjU6gNBEXketEZIuIPBJYf7aIPC4ia0XkC771S0XkSfe1vyjLLqM9+MU5LiaddX1a8u5vwm80jTLDKdcDVwE3eitEZAFwKnC4qg6JyBvd9YfiFIF5E/CHwE9EZI6qvl6ifWOSYApflQRj38Ep7O0QzLTXwcTbaCqlibiq3i0iMwOrPwZcpqpD7jZb3PWnAre4658SkSdxSjD+siz7xjph+dhVCHvQIw97HoV/mzx54nW4kRlG2YiqlndwR8TvUNXD3OUHge8DJ+FU6DpPVe8XkauAe1T1m+521wLLVfXWkGMuxq3s1dPTc1RfX19p9hfB4OAgEyZMaOs516/ZFPna7MNHd9Bpl41Bu8Jsidp2+kFTMtno3z/LeZKIO1YVn3UW6m4flGvjggULVkc1aUjLYYd3aYamEC2fLw3tzk7ZHdgHeCvwFqBPRA4CJGTb0LuL295oGcDcuXO1kzuVROGdLtzbvRcY6aW2y8bLz0gf6gm+h/Nvfu8IG6M8+fAwzb2FDKAm2V33zjl1tw+aYWPdaLeIDwC3qeP+3yciO4Gp7vr9fdv1kqJBqBFP3Uq8thLeWL9mE2m+22UNlBpGXWl3iuG/43SrQETmAF3AVuB24HQR6RaRA4HZwH1ttq2jaKKoNdFmwyibMlMMb8YZmJwrIgMichZwHXCQm3Z4C7BIHdYCfcCjwA+BJZaZUj5VpfFlYfnmq1OnK5aNlaI16kiZ2SlnRLz0wYjtLwEuKcseI5yF05fkyvzwxKxdKY3OgOK9pZ7DMJqITbs3WL9mE5efEZ+CmOSBtiOFMW6qflFsP6Q3dH3XuoHh51Wla1ZN3GCykQ8RGQesAn6jqieLyBeBdwHbgQ3Amar6UtwxTMTHOO0QPo929aMsq7pglIB7rwWFPI1NhpHAOcBjwN7u8gpgqaruEJHLgaXA+XEHsNopHUKSyKR5PRifjqOK2LHftqCdRYhsnMgbRlZEpBf4S+Aab52q/lhVd7iL9+Bk6sViIt5BtBIuSXucvMcsCv+NJsuNp1WBtgFPIwf/AnwG2Bnx+oeB5UkHsXCKMUwZg5R1ruVSNGnGFppK3HWr02dcNi++vie3vnJkyq3vnCoiq3wrlrmTFRGRk4EtqrpaROYH9xSRC4EdwE1JZzFP3BiFP5WuqA49dfFUO0FojNqwVVXn+R7LfK8dC5wiIhtx0q2PFxGv7Mgi4GTgA5qiLoqJeIdRVFw7jyh3ipDX5X0Wid38ikVVl6pqr6rOxKng+lNV/aCInIQzkHmKqv4uzbEsnNKhZA0tFNmjsg4E33/XuoFUcXFvG3+mShhVhBg6NfVxjHEV0A2sEBFwCgN+NG4HE3FjBO1MSWwywZTDutBJ8emxgqr2A/3u84Oz7m/hlA4mGFqJS9EbixT9HmcfPiO0TEAVFH1eu6nXF/PEjdKmtHfyz3t/uKYq7zjq+psgjy3MEzeA4ibDBB8nLLi0AOvaQ9bwiD+GHnb96nADCwp2mQJuN4dqME/ciCTL4GfcoOAJCy5l5V1LizKrUMrOHW+XkKfJJDpt6TFttyMutFSHm9xYwETciKWo+K7nkddRzP1CnjZLxaPqAc66e79p01RN0PNj4RQjNUV80eoWXrEa4fXAPoP8mIgbqSnqi1Y3IQ/SSmzcyM/C6UtYv2bT8I3VhD0dpYVTROQ6nKmjW7xu977XzgO+CPSo6lZ33VLgLOB14BOq+qOybDOyk/SFyhqG8At5HUMsWd9PVYz1vP6kzJ6073ssh2vKjIlfjzP76Eb/ShHZHzgR2ORbdyjO1NM3AX8I/ERE5liLtvqQZgAwr/BVKeje+/LbXcdJPB5+URurwp2WTn//HmW2Z7tbRGaGvPRlnPKL3/etOxW4RVWHgKdE5EngaJwenUZNSOP1eQKYR8xfntXNvI9cEfraWcftx3kfuYJJG4ZGrC9C9JdvvnrEjSTvjeifLvqTlm1JS6cJWCuDoGPZCweQFEWy8h/cEfE7vHCKiJwCnKCq57jVu+ap6lYRuQqnRoBXxetaYLmq3hpyzMXAYoCenp6j+vr6SrO/CAYHB5kwYULVZsSS18b1azYlbwTo+K5U273eLZGvTZ3QxdbB7aPWjxtS5syZlur4cTzxxLMtH2PKlG5eeGGoEHviSHvdg+wzbS9efPa1gq0plqJtdCayOSxYsGC1qs5r5Xj7vWmy/l3fn6ba9h8Pu7Pl86WhbSmGIrIncCHw52Evh6wLvbu45RyXAcydO1fnz59flIml0N/fz1i10dsli1eY5OW+PKs7dP1Zx+3HtT/7TfhO9z3Dqms+ldqGIFkHWj0bg78K3n/GgXz75qeAp0asLzpENH9+Pk/8tKXHcNul9W42XaSNY90D92hnnvgs4EDgIbc6Vy/wgIgcDQwA+/u27QWeaaNtRgtkic/GxZu3H9I7Shg9xh0d/4vRC8PEhVv8uepFZMj4bzj+8748q3vEchmTnfLGxJNuonUeDzDCaZuIq+rDwBu95UA45XbgWyJyBc7A5mzgvnbZZrROVlEJi69HCYgnPEGhDCO4fkSM3RXXMlIcX57VzevdMnz+KG+9SLJe8zRhraonLxnZKTPF8GZgPjBVRAaAi1T12rBtVXWtiPQBj+K0JFpimSnNI226W1i1xKiBq4XTl9C1bgDZtt+I7BdPHKPEPIqgl5yHSRuGYs/76gFOdHDi0+WNN2Ul6yzUIHUR9iw3rk6ZEVpmdsoZCa/PDCxfAlxSlj1G+4iaqh/3RYp6zVv/9a84marBNMYkQS2LpJuIJ+De61XVjykq1z2Lh150umbwf8PSK0ditVOMUvG+cK16Qv5yucE0xjSedVlCHybmYR74y7O6Q8M4ZQp70ZOV0nY1KpK4m7sJuYOJuFE6Rf2UDWupBunEqszYdFqCNxIvPh8caO1aN5D6mpUxYzPMTj9xYl70jSOPAzDWZ7EGMRE3GkWWAdEgdZtG7w+1+Nl+SO8I8UozxpBVuPzXzH9dwn6xRI0jhIm5P9xVlMceJeTBm3pY2KUTMBE3GkkebytLWYAkbzQOb2AzDn/IJSyTZURDjYIGGmXb9uFr4I9xp7kufhuDqZXBGHgZ4ZawzznYCq9piMh44G6cxsi7A7eq6kXua2cDH8dJ9LhTVT8TdRwTcaPRFBkbjYubZ4mpT3xaE4X81QMkMnYeRVDkPUaEZGIGFXX8gSP2CxPuoEgHM22C9oUJelkx87r9kiqAIeB4VR0UkT2An4vIcuANOKVIDlfVIRF5Y9xBTMRrQKekQtWZKjJcgkLeig3zPnIFxHjIee3zP09KmRz21ilOyJvaMSoN6tQ8GXQX93AfCnwMuMytJYWqbok7jol4BFUJayc3Fx4rpAmn+LfNkk8e9HzjJj0FJ0f5JyN5+2ch/fvqTiXkRXjWJyy4NNNAcKu8smM8K549JOXWd04VkVW+FcvcsiHDiMg4YDVwMHC1qt4rInOA40TkEmAbcJ6q3h91FhPxFLS7N6AJ+dghKNBZBB5GpjDGhTK8Y2fx7JMmLA3NGF1wzKN7U/Tsz1cPEF49YDy9jB7gTBJuHd8VuU07ZsEWzNakAljupMYjRGQy8D0ROQxHl/cB3gq8BegTkYM0olqhiXgNsJzXsUWcZx0WL4/yxvPOSo1iZ5ckirPHzN7nOXHauhHr/B7oRnoSjzFw/HgAJs6aFSq8Ye/L+7UQ3N6/bQPFPBFVfUlE+oGTcGpJ3eaK9n0ishOYCjwftq+JeA6Suo3kwYQ8P1kLcNVxgMwv7FGDiFmP4+GJ9s6unUx/W0QlyABBAfevW/HsIczsfZ6NA9FCPrP3efAu89ucdl2bf7Ff6tBRmvII3jZ1/DzTICI9wO9dAX8D8A7gcpw4+fFAvxta6QK2Rh3HRLwFio6b+9PmLKSSjVZvgmkEoyg8IYsKraSJk6cNy/i97u6uHcPP/SIdF+N9z94PcOsrR6Y6V5ATp63jPXs/MLx867QjWfHsIWwc6BkVjonK1Ekr5A1lOnCDGxffDehT1TtEpAu4TkQeAbYDi6JCKdBhIp5FdLOKgie6wX1OW3oMNS8nbgQIE4YyarTEpSJGCXnWmHoUSYNz3uth24V54N2buhLDNMOCPm0dK3pHirn3vryQT5ZfI15HqFZqyleBqq4B3hyyfjvwwbTH6QgRDxPjNIOVeYQ8aX0a7zrtec1bH0lRIalWPLsosc8jvkUJtsfQ9t3Z+IIjwDN7Q8OrQLhIx+EJsfd3aMZ25xghIRm/Zz+z9/kRsXW/dx723pN+ncwLad/XzsyVqugIEU9DU+PRJuQjKbOhcxpvvCgBjyLo7cZliQRf3+3g3YaXkwYm/fslhX/C9h2asZ2vrzoO5jEipBLEH1sfmrGdnV07I7dNU+I3+CvKm/na5HzyJEzE20wW0S1jALUTSNvQOa+Qe1QxQagoujd1hYrzLoEcLZR+8QwK+mhhHSnkQU6ctm44VOP9KvCLuWdjGEli3vA4eWY6QsTrlvmRJM5BW4PhmLDXTehHU/bnHueZFx0KCSPJC4eRQrfbDB2x7AndpA3pzpdUbjfMPk/IPaEOy3rxmNn7PN2/2294eWjG9sR8dBPycjv7XAecDGzxdbv/IvAunBHXDcCZqvqS+9pS4CycbKRPqOqPirSnjuUpg+VH05Am7g7mvXvECXnR6YbeoFwRRHm9YVkdceeMEjr/IKD/+sRdD3/eetocdk/IPS97he+1jQM9sbF5SO+VQ/SNZSyHUqBcT/x64CrgRt+6FcBSVd0hIpcDS4HzReRQ4HTgTTg9Nn8iInPKbNEWFDm/dxwniFHV1MK867gbRvDLYr0NyyPuBp6lJnmQl2d1jxCRnV07R6b0pfCU0xAf9sgu5EEPNeyGH7xWwW5K/udhQu63yZ+5Ehw0DS4PTQiXpCQxj2KsCziU257tbhGZGVj3Y9/iPcB73OenAre4BV+eEpEngaOBX6Y9X5aQQlQmStzrSa9FHdNrKxakFfFIgxXVGk1RXrknWv7Zj04oICBIIeGAVqbh5+3Z6e33Rwfsy6pF70+1T1if0zC2H9IbObvSf/PxX4c0M0ajyCvmYxmJySF3NnBq3p4MHIfjJf838AhOjdu1CfvOBO7wwimB134AfFtVvykiVwH3qOo33deuBZar6q0h+y0GFgP09PQc1dfXx/o1m0Zs47TyGol/m7DXy2JwcJAJEyaMsqGdJL1fv411pUgboz6HNN3gwZkavrNLhjMpvEk0U3buyQu7/Q5w0vkAdtu+2/B+u20f/V3b2RUt4t7244Z0+LxRJB3njw7Yt9TPef2aTaHXz29z0MawTJR9x3Xz3OvOTcE/OcmPd209oq7xuCFlzpxpw8sLFixYnVTLJImJc6fpUf+aLoX7P9/xzy2fLw2xnriIXIwTw+7HaXC4BRgPzAEucwX+027SempE5EKcYuc3eatCNgu9u7hVwJYBzJ07V+fPn8/lZwS9q3tjz79884eymNsS/f39zHdn+3iTftodlw++3+D5z7/5vcM21hX/dWyVuM8hjTf+8qxuth73e9g+Mt/6f/7uCL6154MAw/nY/iyQKKK88bgQiGdH0nG8Y6xa9P5Cr2GQ+fPTXc9gCApGeuafnnAA/zz4NBCTy76n88cfivE8c/81a9rkn7wkhVPuV9WLI167wi1WnsmtFZFFOJ79Cb6ppAPA/r7NeoFn0h4zTXy7TgObWWwponN4Uquv9Ws22axSl7iwSljs1y8kQxN2TaaBdAIO6ZpIhOGPR+eZHFMWUTH2XaGX7ujB2kPTn8efmujdCCY+vUfHiLdHrIir6p0Jr2/B8c5TISInAecDb1fV3/leuh34lohcgROymQ3cl/a4QeLi23WKDyelwIUNfkI+MU+6YXRiDD1vCqIX481bqrUo4gTcT7tELWksKU7IPfyhkbT4Jwy1I7WzbqQa2BSRecCFwAHuPoLTmOLwmH1uBuYDU0VkALgIJxulG1ghIuDEwT+qqmtFpA94FCfMsqTMzJSmU3Qj2iCdJOhhQh7mjYcJZlDId9u+W+hMR484gUnbNSdImqyUutQV8Qs5jM66mfi0stsMjb1BRpXH9YS8lUHTppI2O+Um4O+Bh4HoebE+VPWMkNXXxmx/CXBJSnvGDEkiGeyE7qcdaYmdMJEoq0fuz5AYMa09MJkG0nuGWTv8ZGXeR67grOP2Y35pZ0iH9790woJLR2WxeDjLXWwkOY88S3ncsUra3y7Pq+rtqvqUqj7tPUq1zACcPNe4XNd21FL2SuPWYTwhiN+2Vuxr5UY18enR4g3t+WmfVfjnfeSKkizJxsq7ljJpw9Dww8+kDUNMfNrxyNOKsifmM3ufZ+YNlxdub51J64lfJCLXACtxOjQDoKq3lWKVMYq4Abd2ThSqk2eeVJ2yKDuDoZQi4t3+wcxWPfA8Qu4JZ5WTYVbetXT4l6Znz7ij1bfcTVqPHHxeOTDzhsvp3tTF4//0yTJMrxVpPfEzgSNwWge9y32cXJJNRgjLN1/dUTM6k7zrtKV6WyVp8DCLgHoeu7dPlAdfNE7Py3D748J17SDuJjJpwxC9P902wiNP06T4xGnr+Nt5P2P6237D/JXnMX/leVyyduzKVVpP/I9V9X+Uaolh+PDHqRdOX8L5N783dLukFMykXw5xuc2tpu9lFeg8jSeSqiqOLAsgI7b19vULeR2nqff+dBsvz+pm8wH7pW4v10mk9cTvceubGBUS5403tc9gHFnLHxR5DbxJKVuP+/2ojIeiveewuHDa/fz/D1mPESb6Jyy4tO3eedKNo2vdwHCcfPMv9hvljb9n7wdCHydOWxdbNXGskFbE/xR4UEQeF5E1IvKwiGSapdlplDUQWId4dDsHOL33GzZVvgwh9/aP6grfioD7B/KiBvWyEvd+g78gwqb9R9FuMfcG8GVbeIqgX8g3DvSkCqt4zSjqKuQiMl5E7hORh0RkrYh8zl0/RURWiMh69+8+ccdJG045qWWLO5SiZ0MunL4EauB1V9GwIuycRZST9e8fNi28Fdpd09oLyaR9D0l1tz0hb2eYJerXZte6ASbRy6sHjGcjPazAEehbXzkysntQK42e28AQcLyqDorIHsDPRWQ5cBqwUlUvE5ELgAtwJkmGEuuJi8gEAH9aYTDF0NvGiCaYBueP9RbRw3OsE1VWIew1cEQ5rbD7499+8SvCC6+iKUFUWdgo29PaWPUAqEfXuoER6YeeR15joY5EHQbdxT3ch+JUdb3BXX8D8O644yR54t8XkQeB7wOrVfU1ABE5CFgAvA/4OjCq2qART6cKcl5mHz6D5Zs/NOIGOExM6iUke5NhnmtSCVkYLYDjjta2Cne7u9fUrVdl9yYn/fDrgeYScT09W2X79t2zTCiaKiKrfMvL3AJ+w4jIOGA1cDBwtareKyL7qupmAFXd7NaoiiSpdsoJIvJO4P8FjnVjMzuAx4E7gUWq+mzad2SEN42oU+513Wml5Vpw+rl/1qBHWCPi3p9uy3W+POTJUMm7/1hoX+bvHBQU8xqwNakUrVte5AgRmQx8T0RGle1OIjEmrqr/AfxH1gMbI4lqGmEeeWvEhU3CxGzEjMUQAfeLgN/j8jzQdoQVsghxlDeeFPppqoD7664EJ18FxbyuA5phqOpLItKPM/74nIhMd73w6SQUGcxeMsxIRRbP2rzwbERdLy+u7T08vMkucY+kwknezXblXUtHDLx1rRsYfsyZM204y6LVsEOedEH/+y5LpMu+ibXSsMWrZeOPl6fJYqkKEelxPXBE5A3AO4B1OFVdF7mbLcIJZ0diIl4iacS5qQK+fPPVldsezCrx8ItzGvwC7nnf3t8wj3b55qsLEeok0gpxmNce58m3KvBl/3pM+3+VNOM1yyzPipgO3OWma98PrFDVO4DLgBNFZD1worscSZmNko0C8X9xotLqyi5RG7TH32iiCkEPuw5B4U4Scn8IxS/g3k/0SRvi4+H+993f3z/iNX9tkCrx1yVpShglbOwoKnQWJuRD7evAmBu3I9qbQ9b/Fjgh7XHS1hOfBQyo6pCIzAcOB25U1ZfSnqjT8P7xTlt6zIjlrGKXx+uJixMXKfB+29ot5GHjCVGCHdaBPq6bjteRxxO8VjoqtZrT3+pAJ8R75Vmpon6P91l7N+00reqSmnaMJdJ64t8F5onIwTg1wW8HvgW8syzDmkgrRZnCPI84vC9T1okuZU7Pb3czieWbrx7OMMk6QSdNHNx/nrwUMRmpCCEfK/ivp/874Am7N9D76gHjK7Ox3aSNie9U1R3AXwH/oqqfxInnRCIi14nIFhF5xLcucjqpiCwVkSfdqf1/kefNVE0rX/a8NbGL8Ix0fFemCTJpaEfWTViKoEdcrDT42saBnuFQStq+mGko8kZWhzBIlVU0/WMw3kCyR1j9mKr6i1ZBWk/89yJyBs5I6bvcdXsk7HM9cBVwo2/dBYRMJ3WLa50OvAmnx+ZPRGRO01q05RWuVlMNi/D2PIqsTV7m1Hy/Bx41wLVr3ciJO3FhlKLxd7JpFb+3mWc/jm7ZhMoJ+8Ua9v1Zddfo/7l/ZGz230wr4mcCHwUuUdWnRORA4JtxO6jq3SIyM7D6VBjuEHUD0I9TE+BU4BZVHQKeEpEncf7lfpnSvkYRlTPukUfQi6xuWPQAaVmNGjwxDrY2Syt2YXFT7zhFe75FDnK22yuvax37uIbonYSopvvZ4eYxzlDVx1Mf3BHxO1T1MHf5JVWd7Hv9RVXdR0Suwmma/E13/bXAclUdNZ1fRBYDiwF6enqO6uvrS2tOW/Gq7u0zbS9efPa14fVp82DDqva1go6P9jKnTOnmhReihSGqslxe8uQCDw4OMmHCrjI9azePnP/gVegbNzTy//n1bhlRRzsOf5W/4HEA5syZlsnGMJ54oroJzkmfc5C4z72VfO440lzDvCxYsGB10gzKJLoP7NXpF5+datun/+aCls+XhrTZKe8CvoTzu/RAETkC+F+qekpBdoR9y0LvLm7tgWUAc+fO1flFlggsEM+sr3/lRm679N7MnsLlZxQfU47yyt9/xoF8++anYvct1hvLfj36+/vxPut5H7liVEhk4tNKz50bRu23/ZBeBo4fPciVphZKkK51PwOivT6/jVF8/nPVpRym+ZyDRH3uyzd/qAiTRpHmGhojSTuweTFOeOMlAFV9EDgwx/mec6eREphOOgDs79uuF3gmx/Frh1O4qR4/9fyzC7NS54FPfzpg2u1hV33vnjs3pBDwYm5idQ1NRBH2mdfl/9lwSCviO1T15cC6PMO/UdNJbwdOF5FuN94+G7gvx/GNlESJSdi0dT9Fi3lWwrzwOAHuWjcwqoCVt7/3Hr3B3LAbnDcz04RrF1bvp16kHdh8RET+JzBORGYDnwB+EbeDiNyMM4g5VUQGgItwpo/2ichZwCbgvQCqulZE+oBHcaokLmlaZkoT2SVY2X9UJfW2TEPWCUJRYZQkG7rWDTBx1qzQ7vJeXrE/Kyc6hFCMkBeZTWQYaUX8bOBCnE4U3wJ+BHw+bgdVPSPipdDppKp6CXBJSnvGNO32dIoeuMxCFiEPE/Bgn8koHG89/NfFcNEoettWL7vpQm7lk+tD2nDKX6rqhar6Fvfxj0BRg5qGj6p+qnpCmCd9rVUxSvOeg9koHml/BXg9Gv0EKxmGDYAa0ZTVR9bIRloRD3NP6tPiY4xQ9RfCL+RZxbzVWHncex9RA5xd1euy2ugX8qgJP3M//+VMx8xKsD9o3WmCjZ1ObDhFRBbi1EfZT0S+4ntpb5zYtVESrU7+yUvWn/nBAdBJ5BfyyMkwgXO8eoDQ+9NtuQTG32w3Ck/IH/+nT2Y+fhr8MwybHlYBC61UTVJM/BlgFU7oZLVv/atAOf/hHUpQqP1fjHZ3AEorLMEa3s7f8aXUrfBPxCkCz8a4wllzP//lEfnnZQlVHYS8VY/bhLw6knpsPgQ8JCLfUtXft8mmjiTJ8/bX7vavK0vck4QlTMCjlmF06deoUrBpbgDeNnnLA3StG4BZs0adL6kSYpFCFfw881alNIy0MfGZInKriDwqIr/2HqVa1sHECYXfO/f+BmtIlN11J09Z1LTNGuKENKrQVdZ4fKux+yJvnMHPyWLQRlbSphj+G06e95eBBTgFscZmSbCakEbI02zfSmy9Xd5h1vBL1IBmq7nrUb8Owio7FlnUK/iLqt1eud04mk1aEX+Dqq4UEVHVp4GLReRnOMJu1JSguISFZNKQFFqJii97FQL93XTCSJNlkrW1WNDeYP3pSfSO+kURDPektc3fxSlv2Y+oG0HcZ2Whl2YjIvvjlOqeBuwElqnqlb7XzwO+CPSo6tao46QV8W0ishuwXkQ+DvwGeGNe441mE1XqNUzMgzW6k6bIB+Pcnrh9/Ss3jvIYk0RshH1uDNx/fq9jTtQsUL9daSl6gC/uWAunLzEhbzY7gE+r6gMiMhFYLSIrVPVRV+BPxJnZHktaET8X2BNnuv3ngePZVQPFaBh5BkTDfuIHBdkTTS8s4RfwpLxu/5T3MOFyColFV85LW6s7bTy/1f6a7cD7HE3Im4mqbgY2u89fFZHHgP1wyo98GfgMu+pLRZJqYFNV71fVQVUdUNUzVfU0Vb0nv/lG1eT1FuMqIXqThCZtGBoehAwORvr3DztWXruC0+XjxNo/SzOMQhokt4lWvf6yGmcbgFM3apXvsThqQ7f3wpuBe0XkFOA3bnZgImnric8B/h44wL+Pqh6fZn9jbBI3AFdFT0i/kPs9VG99VE9O75dD8GbTCu1uGp2XrF583HXZfkgvJyy4NPLX1FhAtkuWNn5b0zSFEJEJOM3oz8UJsVwI/Hnak6RNMfwO8ADwjzhi7j2MBhOVmpj2C1h1Wdo4lm++epSX74l52C8EL9yTtqBW2PmiaseX7aG2q2gXpPvMtx/Sa155SkRkDxwBv0lVbwNm4ZQVfUhENuL0VnhARCLbSmWpJ/5VVb1PVVd7jxbtN2qAJz5h6+PwvshJ9cej9oui6Bzs4PtYdc2nWHXNp0aEfoK1YrLemKKygNpJO4Uc0l0jE/J4RESAa4HHVPUKAFV9WFXfqKozVXUmTsOcI1U1sq9fWhH/gYj8nYhMF5Ep3qPVN2GMLdKIeRHx7yJYedfSQuLBafLzoT2CVrWQh11PE/JYjgX+GjheRB50H+/MepC02SleJoo/hKLAQVlPaBjBME67v+je+VoJBaW5AeXNy2+FlXctTZ2pUwZhcwqsrko4qvpzEiZNut54LKlEXFXz9NOMREQ+CXwE50bwMM4M0D2BbwMzgY3A+1T1xSLPa6QnSXiGy9aGTJrxuuWE0YQvs39gromepH8gt0yifsnEeeRN+PybRmw4RUSOd/+eFvbIc0IR2Q8n33yeqh4GjANOBy4AVqrqbGClu2xUQBHCFRZWqUvO9fLNV4d64WEhljwzKaNo9w3B6w9aFnl+yVgjieJJiom/3f37rpDHyS2cd3fgDSKyO44H/gxwKnCD+/oNwLtbOL7RAlmyVMI65kTRrkyWNCIRFDdvOex9R12LJgg5EJnXH7ZNXW60RnpEtfjaz4knFTkHp5/mfwM/VtUPiMhLqjrZt82LqrpPyL6LgcUAPT09R/X19bXJ6nwMDg4yYcKEqs2IJcnG9WuiZ/7qeCdn9vXu8NDeuKGR/19z5ozOlAoePyxbJst19I4Xdhw/Tzyxa8Bftm1P3N5/bD/eflE2xu3TLtav2cQ+0/bixWdfC7UlaKP3uSbRan/W4HUo8/uyYMGC1WnytuMYv9/+esBHP5Vq2yc++6mWz5eGWBEXkVhrvbSYTCcU2QcnL/L9wEs4Oei3AlelEXE/c+fO1ccffzyrCW2lv7+f+XmrIrWJNDZGeZB+7zoshBL00sN+3gePHeb1ZrmOSfHXsKnqaSeoxNV6j7Ix6tq1Oz789a/cyG2X3ht5/jA7k349teq5B20o8/siImNSxJPCKRPdxzzgYzjz+vcDPgocmvOc7wCeUtXn3UYTtwFvA54TkekA7t/wzrjGmKMsMQsTpbB1eYWo7LrtReP3etPaHXdtigi9WHy8dZI6+3wOQER+jJNw/qq7fDGOB52HTcBbRWRPnHDKCTgt4F7DSWW8zP2bWPjFaB9pUgG9qoDec4+wdWHHL4M0WTZlx+qrSKOMoo43HUtBbI20k31mAP7g13acVMDMqOq9OOGTB3DSC3cDluGI94kish6nBONleY5vlEeaL1pw9qOfPB2B2oE3oJdVaLNsHzZYXBdh94j6fIMetw2A1ou0Iv4N4D4RuVhELgLuZVcmSWZU9SJVPURVD1PVv1bVIVX9raqeoKqz3b8v5D2+0T6a+mXOK0StzsQMa4hdJzGPE3IT73qSKOLu/P4bcSbkvIgzGHmmqlY3LcwwSiKPqKbdp05iHUUVNloopTUSZ2yqqorIv6vqUTghEMMYQVxcua4hlFbIM5szLqulDlR1g6nTNWgqaWun3CMib1HV+0u1xqg9VRR2qivBAcv1azZF9tiss1hFpXh28mfbJNLGxBfgCPkGEVkjIg+LyJoyDTOaQZw4hXnhVRZnKoOm3tS8EFCweUVYjXmj3qQV8YU4FQuPZ9eU+3eVZZQxNqmi209essTGmyZ0WW80TXt/nUbaHptPA5PZVTdlsrvOMMY0WYQ8qrNPXVi/ZlPuXwpV5fEbyaQScbfWyU3AG93HN0Xk7DINM8YeXnjFvrjtJexXRTBvvc43HyOetAObZwHHqOprACJyOfBL4P8ryzCjOSzffHXmWLfN0iuXuBtl3nZyRc88tc+/GNLGxAV43bf8OgkdKYzOImwSSJNi4GMVv8ddRLinSOG1X2TFkNYT/zfgXhH5Ho54n4rT4NMwgGze+PZDekud+Vfnqe3twrzc+iMi1+EkiWxxG+QgIkcAXwPGAzuAv1PV++KOk3Zg8wqcGZsvAL/FmbH5L3mNN8Ym7W7Um4Yi0uQ69UYQRR1SD/0pkg3+fK4HTgqs+wLwOVU9AvisuxxLWk8cnBCKuo+dGfYzjFG0q8tPUWTtERk1Q7PuMzebQJr6801AVe8WkZnB1cDe7vNJOF3PYsmanTIVy04xMjDW4uJ5a6TEZYE02JMslIXTl0R2kUozOamGTBWRVb7H4hT7nAt8UUT+C/gSkPjz1rJTjEJpR33uOpNGkMNi9p2arZO2cUddGLcdJj6duqXl1hydfT4GfFJVvysi78MZe3xH3A5pRdyyU4zUpBXydgpXkelxYXYvnL6E05YeE3reJLuMdHTItVoEnOM+/w5wTdIOebJTwOlEb9kpRmZentU9IsSSNdZcF5JuCE17P3Wng67nM8DbgX6cMifrk3ZIJeKqeoWI9AN/iuOBn6mqv8prpYhMxrnDHIYTyP8w8DjwbZyOQRuB96nqi3nPYRhByq7O5+Rhf6iUY3cqY1m8ReRmYD5O7HwAuAj4W+BKEdkd2AYkxtFTibiIvBVYq6oPuMsTReQYt9VaHq4Efqiq7xGRLmBP4B+Alap6mYhcAFwAnJ/z+EbF1DU2XpaAL998Nf39/aUcu66UfVMcywIOoKpnRLx0VJbjpJ2x+VVg0Lf8mrsuMyKyN/BnuOEYVd2uqi/hTCDyWr7dgBOyMRpKJwh4WN/MTiTv+4/br9OvaRZENXmkVUQedJPP/evWqOrhmU/ozEhaBjwK/DGwGieQ/xtVnezb7kVV3Sdk/8W4PzF6enqO6uvry2pCWxkcHGTChAlVmxFLkTauX7MJHd81vPx698jx73FDzv+bbNs+Yv3sw2eUbmNU+loewuyt+2ddtn1Zrq//+vn3m37QlNJsXLBgweoc2SIj2Ktnf/2jUz+ZatvV13665fOlIe3A5q9F5BPs8r7/Dvh1C+c8EjhbVe8VkStxQiepUNVlODcB5s6dq/OjWqnUhP7+fjrFRs/T9XvhwcYQ3qDm6Gn398Z6X0XYePkZRXrio+2t+2ddtn1Zrq9/7MDbzwtJ1fka1pG04ZSPAm8DfgMMAMeQIuAewQAw4Iun34oj6s+JyHQA9++WnMc3KiBMwGGXaE/aMJQ48adpU6ibZGudCKuiaOGT/KTNTtkCnF7ECVX1WRH5LxGZq6qPAyfghFYexcmRvMz9+/0izmdUS9ATL7PwVTuw/pNG3chSOwUAEXlAVY9s8bxnAze5mSm/ximutRvQJyJnAZuA97Z4DsMojLw1uMci/tx+u5lVT6yIi8h/4JRC3Ohf3epJVfVBICzgf0KrxzbqQ3BiTxztEMW8otPJgp2ECXn1JMXErwd+LCIXisge7ro7yzXJGEuEdbw3mk0whm03uWqJFXFV7QPejFMacZWInAe8ICKfEpFPtcNAo9n4PfGq4+FppsqbIBlNI012yu9xJvd0AxMDD8OIpK5laNMUpfILuoULksly82taFlLdSYqJnwRcAdwOHKmqv2uLVUbtCVby82KjcdPto9qytbOSYd5zdmqp2LTkEWW7psWQ5IlfCLxXVS8wATeCJH1xX57VPSomHhT4un+J627fWMA889aI9cRV9bh2GWI0kzRfvqgslaYIZFPsrJKsWSo2MFocaWdsGsYwSV9WL2TizdKMio2b99WZmGgXS+bJPkZnE+xxGLY+iD+kUtfBTqM9mIAXj4m4kYukL2Nd64kb5RB2c4+64RvFYuEUIxNRudT2Je1MwjrQR21nlIOJuFEYeYTcvtzNwxPu4GdnN/JqMBE3SsWm3XcGUQJus2CjEZHrRGSLiDziW/dFEVknImtE5HtuP+JYTMSNUvEPZFqWSvMJC53ECXjYc2OY64GTAutWAIe5XdOeAJYmHcRE3Cgdf5phlTM2jdbIIuBGMqp6N/BCYN2PVXWHu3gPkJgdYCJuFEZej9o88WZggp2ZqSKyyvfI2g3tw8DypI0sxdAoBBPizqDThXzcNs0y12Fr3kbJInIhsAO4KWnbyjxxERknIr8SkTvc5SkiskJE1rt/R3W6N5pPWO54pwtD07BaJ+UiIouAk4EPqKombV9lOOUc4DHf8gXASlWdDax0lw3DqCkm5sXjVo49HzglbdHBSkRcRHqBvwSu8a0+FbjBfX4D8O42m2W0iaA37s87NlGoP/bLqRhE5Gbgl8BcERlw+wtfhdOrYYWIPCgiX0s8TgpvvXBE5FbgUhxjz1PVk0XkJVWd7NvmRVUdFVJxBwcWA/T09BzV19fXJqvzMTg4yIQJE6o2I5YibFy/ZtOodTq+K3J72bY98rXZh88Yta5TrmOZFG2f95mHfV55KfMaLliwYHXeGLXH3hN79S3z0jkaP+3/h5bPl4a2D2yKyMnAFlVdLSLzs+6vqsuAZQBz587V+fMzH6Kt9Pf30wk2Xn5G+D920OtOatEW5eV1ynUsk6LsGz1T80MtH9Oj7tewjlSRnXIscIqIvBMYD+wtIt8EnhOR6aq6WUSmA1sqsM0omKr7ahrFYaGuetL2mLiqLlXVXlWdCZwO/FRVP4jTAm6Ru9ki4Pvtts2oFou1GkZ26jTZ5zLgRBFZD5zoLhsNoVUBNgGvN+aF15dKJ/uoaj/Q7z7/LXBClfYYhjGSOPG2G289qJMnbjSYVj018/TqR1ipWauZUj9s2r1RGzzRMIGonrhGD/b51AvzxA3DGIG1VWsWJuKGYYRiAt4MLJxiGMYITLybhXniRq0wATGMbJiIG4ZhNBgTcaNlLD3QMKrDRNyoFWElaRdOX8L6NZvsZmEYIdjAppGbskXVRNswkjFP3MiFCaxh1AMTcSMzJuCGUR8snGJkomwBj2sKUWTzAcMYK5gnbhiGUREiMllEbhWRdSLymIj8SdZjmCduGIZRHVcCP1TV94hIF7Bn1gOYiBuGYVSAiOwN/BnwNwCquh2I7iAeQRWNkvcHbgSmATuBZap6pYhMAb4NzAQ2Au9T1RfbbZ8RTd4GAWni6Dbd3uhADgKeB/5NRP4YWA2co6qvZTlIFZ74DuDTqvqAiEwEVovICpy70UpVvUxELgAuAM6vwD7DJe0gpgmw0SnItu1Zmn9PFZFVvuVlqrrMt7w7cCRwtqreKyJX4ujeP2Wxqe0irqqbgc3u81dF5DFgP+BUYL672Q04bdtMxGuOCbhhRLJVVefFvD4ADKjqve7yrTginglR1TzGFYKIzATuBg4DNqnqZN9rL6rqPiH7LAYWA/T09BzV19fXHmNzMjg4yIQJE6o2I5YwG9ev2ZRq39mHz0jcJs2xko7T1OtYJ+puH5Rr44IFC1YniGoik/Z4o75t6ntTbfvDZ/818Xwi8jPgI6r6uIhcDOylqn+fxabKBjZFZALwXeBcVX1FRFLt5/4cWQYwd+5cnT9/fmk2FkF/fz9NsjFLHnhaL9x7+/Ex9fgc8KZdxzpSd/ugGTYWzNnATW5myq+BM7MeoJI8cRHZA0fAb1LV29zVz4nIdPf16cCWKmwzqsFmgRqdiKo+qKrzVPVwVX13nmSOtou4OC73tcBjqnqF76XbgUXu80XA99ttW6cT1d28XXFvE3LDyE4V4ZRjgb8GHhaRB911/wBcBvSJyFnAJiBd4MkohDABj1tuleDxTMANIx9VZKf8HIgKgJ/QTluM0bTD6w47h2W5GEY+bMamAVQn3oZhtIYVwDKsa45hNBgTccMwjAZjIm4YhtFgTMSN4dmSwQbFhmHUHxNxwzCMBmMiboyobWIZJIbRLEzEO5xg+MTCKYbRLCxP3ADMAzeMpmKeeAfjed1pyskahlFPTMQNm+xjGA3GRLyDsRCKYTQfE3EDMEE3jKZiA5sdQly4xGLihtFcTMQ7GM/77u/vr9YQwzByYyLeAXheuIVMDGPsUbuYuIicJCKPi8iTInJB1fY0Hcs6MYx6UpTW1UrERWQccDWwEDgUOENEDq3WqmZj3rdh1I8ita5u4ZSjgSdV9dcAInILcCrwaKVWNRwTcsOoHYVpXd1EfD/gv3zLA8Ax/g1EZDGw2F0cEpFH2mRbXqYCW6s2IgGzsRjqbmPd7YNybTyg1QO8suP5H/3w2X+dmnLz8SKyyre8TFWXuc8TtS4tdRPxsAbKOmLBuQjLAERklarOa4dheTEbi8FsbJ262wf1t1FVTyroUIlal5ZaxcRx7kb7+5Z7gWcqssUwDKMsCtO6uon4/cBsETlQRLqA04HbK7bJMAyjaArTulqFU1R1h4h8HPgRMA64TlXXxuyyLOa1umA2FoPZ2Dp1tw+aYWPL5NC6SEQ1VxjGMAzDqAF1C6cYhmEYGTARNwzDaDCNFfE6Ts8Xkf1F5C4ReUxE1orIOe76i0XkNyLyoPt4Z4U2bhSRh107VrnrpojIChFZ7/7dp0L75vqu04Mi8oqInFv1NRSR60Rki39eQtx1E5Gl7v/m4yLyFxXa+EURWScia0TkeyIy2V0/U0T+23c9v1ahjZGfbRXXsXGoauMeOAMBG4CDgC7gIeDQGtg1HTjSfT4ReAJnSu3FwHlV2+fatRGYGlj3BeAC9/kFwOVV2+n7nJ/FmaRR6TUE/gw4Engk6bq5n/lDQDdwoPu/Oq4iG/8c2N19frnPxpn+7Sq+jqGfbVXXsWmPpnriw1NWVXU74E1ZrRRV3ayqD7jPXwUew5mZVXdOBW5wn98AvLs6U0ZwArBBVZ+u2hBVvRt4IbA66rqdCtyiqkOq+hTwJM7/bNttVNUfq+oOd/EenHzkyoi4jlFUch2bRlNFPGzKaq3EUkRmAm8G7nVXfdz9SXtdleEKnFlhPxaR1W4JA4B9VXUzODci4I2VWTeS04Gbfct1uYYeUdetrv+fHwaW+5YPFJFfich/ishxVRnlEvbZ1vU61oqminhhU1bLQEQmAN8FzlXVV4CvArOAI4DNwD9XZx3HquqRONXTlojIn1VoSyTuBIhTgO+4q+p0DZOo3f+niFwI7ABucldtBmao6puBTwHfEpG9KzIv6rOt3XWsI00V8dpOzxeRPXAE/CZVvQ1AVZ9T1ddVdSfwdSr8Saiqz7h/twDfc215TkSmA7h/t1Rln4+FwAOq+hzU6xr6iLputfr/FJFFwMnAB9QNNrshit+6z1fjxJvnVGFfzGdbq+tYV5oq4rWcni8iAlwLPKaqV/jWT/dt9ldAJZUXRWQvEZnoPccZ9HoE59otcjdbBHy/CvsCnIEvlFKXaxgg6rrdDpwuIt0iciAwG7ivAvsQkZOA84FTVPV3vvU94tS0RkQOcm38dUU2Rn22tbmOtabqkdW8D+CdONkfG4ALq7bHtelPcX7urQEedB/vBL4BPOyuvx2YXpF9B+GM9j8ErPWuG/AHwEpgvft3SsXXcU/gt8Ak37pKryHODWUz8HscD/GsuOsGXOj+bz4OLKzQxidx4sre/+PX3G3/H/d/4CHgAeBdFdoY+dlWcR2b9rBp94ZhGA2mqeEUwzAMAxNxwzCMRmMibhiG0WBMxA3DMBqMibhhGEaDMRE3EnGrMz4lIlPc5X3c5Za7h+e0Z9D9+4cicmsLxzlXRPYsyKZ3i8hnM+7zk5qUDzAajKUYGqkQkc8AB6vqYhH5/4GNqnppG867u+4q4OStG1TVCQUceyMwT1W3FnCsX+BMqEl9LHcmZa+qXtLq+Y3OxTxxIy1fBt4qIufiTGoKrV0iIh9yCxk9JCLfcNcdICIr3fUrRWRGwvrrReQKEbkLuNydmftLEblfRD7vO9dMry61iPyNiNwmIj8Up773F3zbfVVEVolT4/1z7rpPAH8I3OWeBxH5c/c8D4jId9waOIjIZSLyqGvnl0Le8xxgyBNw1/6vilNb/tci8na3sNNjInK9b9fbcWamGkZ+qp5tZI/mPIC/wJmRemLE62/CmVk31V2e4v79AbDIff5h4N8T1l8P3IFbOxpH7D7kPl8CDLrPZ+LWpQb+Bmfa+CRgPPA0sH/AjnFAP3C4u7zRZ+tU4G5gL3f5fOCzwBT3PXm/WieHvO8zgX/2LV+PUx5ZcMqpvgL8DxynaTVwhG/b9cAfVP3Z2qO5D/PEjSwsxJkyfVjE68cDt6rrkaqqVzf6T4Bvuc+/gePJx60H+I6qvu4+P5ZddVS+EWPfSlV9WVW3AY/iNJMAeJ+IPAD8CudGc2jIvm911/8fEXkQpxbKATgCvA24RkROA34Xsu904PnAuh+oquJMJ39OVR9Wp8DTWpybj8cWnF8EhpGL3as2wGgGInIEcCKO2P1cRG7B+f/5gbvJ13A8zzSDLFHb+Ne/lnIfP0O+568Du7uFk84D3qKqL7rhjPEh+wqwQlVHhTdE5GicBhWnAx/HuVn5+W+cXwBhtuwM2LWTkd+78e7+hpEL88SNRNzqjF/FqY++Cfgi8CVV/S9VPcJ9fA2nCNT7ROQP3P2muIf4BY4AAnwA+HnC+iD/J7BdFvbGuSG8LCL74vya8HgVp40eOF1vjhWRg13b9xSROW5cfJKq/gdwLk7N6yCPAQdntMu7rtNwwjqGkQsTcSMNfwtsUtUV7vK/AoeIyNv9G6nqWuAS4D9F5CHAK8f7CeBMEVkD/DVwTsL6IOfgNLC4n9Eebyyq+hBOGGUtcB3ODcFjGbBcRO5S1edx4uo3u/bcAxyCI/J3uOv+E/hkyGnuBt7sinIWjgLu0UD2jWFkwVIMDaMARORKnDj4TzLuc7uqrizPMmOsY564YRTD/8apg56FR0zAjVYxT9wwDKPBmCduGIbRYEzEDcMwGoyJuGEYRoMxETcMw2gwJuKGYRgN5v8CR6ro2nomSwcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "elev = getvariablevalues(ncdataset, 'Elevation')\n",
    "water_depth = elev + depth\n",
    "plt.contourf( water_depth)\n",
    "plt.grid()\n",
    "plt.xlabel('X-coordinates (m)')\n",
    "plt.ylabel('Y-coordinates (m)')\n",
    "plt.colorbar()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 1.0, 'Vertical section through red line')"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXEAAAEWCAYAAACQdqdGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABDUklEQVR4nO2de9wcZXn3v78E8gRISKB5lEiAcAwFiggRrIeagFiwCL7UA9RDVChVEUWlQqRvpbUU8IDaSvWNQgFBNCIqQqPEyFNqFTBBCGcDksRIICDHgEkIud4/ZiaZZ56Z2Znd2d3ZZ6/v57Of3Tndc809u7+597qv+7plZjiO4zi9yZhuG+A4juM0j4u44zhOD+Mi7jiO08O4iDuO4/QwLuKO4zg9jIu44zhOD+Mi3qNIeqekGyooxyTtVYVNTZ5/gaQ5HTrXcklv6MS5GtHJepc0PTzfVhnbN9eLpE9J+kYn7HKqwUW8DUj6iaR/Tll/nKRHsn5MOeWN+BGa2ZVm9sYq7O0Uks6RdEV8nZkdbWaXteFcl0r6l6rLHe2Y2b+a2cndtsMpjot4e7gUeLckJda/G7jSzDYWLais4DvVUOd6r7NtTudxEW8PPwB2BF4XrZC0A3AMcLmkMZLOkvSgpD9Imi9px3C/qNV9kqSVwM+Am8JinpK0VtKfS3qvpJ/Hyt9f0kJJT0h6VNKnwvWHSvqlpKckrZb0FUnjilxEeI7fSnpW0kOS3hnb9n5J90p6MvznsVueLZKOAj4FvCO8hjvCfYcknRx+HiPpHyStkLRG0uWSJiXqZY6klZIel3R2ht2nAO8EPhme60exzQdJWirpaUnfkTQ+PGaWpFWSzpT0CPCfkgYkfUnSw+HrS5IGYnXz88R5N7tIJP2JpB9JekbSryT9S3J/4A2SloV1eFHKQz8q9xxJV0u6QtIzwHslTZJ0cXhPfx+WPzbcf6ykz4d19Fvgr3Jv9MhzXVGkzvO+x04HMTN/teEFfB34Rmz574Dbw8+nAzcD04AB4P8BV4XbpgMGXA5sB2wTW7dVrLz3Aj8PP08EVgOfAMaHy4eF2w4BXgVsFZZzL3B6rBwD9kqxfzvgGWBGuDwV2D/8/BbgAeBPw3L/AfhFAVvOAa5InGcIODn8/P6w3D2ACcA1wDcT9fL1sE5eDqwH/jSj/i8F/iWxbjlwK/AygofsvcAHwm2zgI3ABeE92Qb45/A+vQQYBH4BfCZZ/2l1CXw7fG0L7Af8Lr5/uO91wGRgV+Ax4KiMazkHeCGs9zGhbT8g+N5sF9p3K/B34f4fAO4Ddgmv80YS35+UenlD8h41qnNyvsf+6qDWdNuA0foCXgs8DWwTLv8v8LHw873AEbF9p4Y/0q1iP5w9YtujdVkifiLw64J2nQ58P7acJ+JPAX8dXUNs2wLgpNjyGOB5YLc8W2gs4ouAD8W2zUipl2mx7bcCJ2Sc61LSRfxdseXPAl8LP88CNgDjY9sfBN4UW/5LYHmy/pN1CYwN7Z4R2/YvjBTx18aW5wNn5dTbTbHllxKI6TaxdScCN4aff0b4cAqX35j8/qTUS56Ip9Y5Od/jbv/++unlvrU2YWY/l/QYcJykW4FXAseHm3cDvi9pU+yQFwl+nBG/K3G6XQgEZwSS9gEuBGYStAq3ApYUsP85Se8AzgAulvS/wCfM7L7Q/i9L+kL8VMDOebYU4GXAitjyitDeeL08Evv8PEGLvQzJ418WW37MzNY1sCe+fxaDBHbH72Ha/SxzLfHjdwO2BlbHPDBjYvu8LLF//BqaIcvOvO/x71s8p1MQ94m3l8uB9xB0aN5gZo+G638HHG1mk2Ov8WYW/+Jbxuc0fgfsmbHtqwR/rfc2s+0J/NKpvtckZvYTMzuSoIV1H8Hf6uh8f5ewfxsz+0UDWxpdx8MEwhCxK4GL49H03fPNr+CYNHseDj8/R/BQBEDSTrH9HiOwe1ps3S5N2JNl2+8IWuJTYvW/vZntH25fnTjfri2eO4si32OnzbiIt5fLgTcAfwvEw+i+BpwbdQZKGpR0XE45jwGbCHzFaVwH7CTp9LAzbqKkw8JtEwl822sl7Qt8sIjhkl4q6VhJ2xEIxlqCVlZk/1xJ+4f7TpL0tgK2PApMl5T1vbsK+Jik3SVNAP4V+I6ViOaJ8SjZ9VWUq4B/CO/PFOAfgShE8g5gf0kHhZ2j50QHmdmLBP78cyRtG9b7e1q0ZTNmthq4AfiCpO3DDsY9Jb0+3GU+8BFJ0xR0qJ9V1bkTlP0eO23ARbyNmNlygs6w7YBrY5u+HC7fIOlZgs6hw0YUsKWc54Fzgf9VEGXyqsT2Z4EjgTcT/PVdBswON58B/A3wLEFL+jsFzR9D0Dn5MPAE8HrgQ+H5vk/QAfjtMFriLuDoArZ8N3z/g6TbUs55CfBNgmich4B1wGkF7U1yMbBfWF8/aLKMfwEWA0uBO4HbwnWY2W8IOj5/SnCNyciTDwOTCOrgmwQPhPVN2pHGe4BxwD3Ak8DVBP+YILjPPyF40NxG8EBpB6W+x057UNgh4ThOG5F0AbCTmXVkdKrTP3hL3HHagKR9JR2ogEOBk4Dvd9suZ/ThIu447WEigRvjOQIf9ReAH3bVIqd2hAOzfi3punD5IEk3S7pd0uKwAZBfhrtTHMdxuoOkjxOE/25vZscoSGr3RTNbIOlNwCfNbFZeGd4SdxzH6QKSphGkRIhnjTRg+/DzJLaEtGbS04N9Jk+ebHvt1bUsqrk899xzbLfddt02I5VO2rbsjpXDlvd+eXbIstdZc9TVtrrZtWTJksfNbLCVMl43a7w9+cSmxjsCd9/5wt0EEVYR88xsXmz5S8AnCVxvEacDP5H0eYJG9qsbnaenRfylL30pixcv7rYZqQwNDTFr1qxum5FKp207euqpmz8vWHxR5n5eZ81RV9vqZpekVkeu8uQTm/je9VMK7bvvrqvXmdnMDFuOAdaY2RJJs2KbPkiQnuN7kt5OECqbmwPf3SlO21mwOlu4HadPeQ1wrKTlBInSDg+zR85hS1z/d4GGHZsu4k5HiIQ83ip3nH7FzOaa2TQzmw6cAPzMzN5F4AOPRt4eTjCQLJeedqc4vcnRU0+tVet8mLunRnY5fcnfEiSX24rAn35KowNcxJ2OsWD1RbVqiSdtcQF3uoGZDRGkZMbMfk4wB0Bh3J3iOI7Tw7iIOx2jTq3wNOpun+Ok4SLudIU6CKa7T5zRgPvEnbYQiXReVEpy3fFzD6OTYcV1eJA4Tqt4S9xpK0dPPbWnxLKXbHUccBF3WiQS6W5HelQpvi7kTi/h7hSnEtJEu8gAn/hxQ0NDTZ27atF1X7nTS3hL3KmEokJdZH27yXPxuIA7vYaLuFMZZYW8iha0uz6cfqdtIi7pEklrJN2VWH+apPsl3S3ps7H1cyU9EG77y3bZ5VRLUpyzRLXs+iI0c6yLvjPaaKdP/FLgK8Dl0QpJs4HjgAPNbL2kl4Tr9yNIArM/8DLgp5L2MbMX22jfqKAOeT+Svu+kTe0WziLX7eLtjFba1hI3s5uAJxKrPwicb2brw33WhOuPA75tZuvN7CHgAQqkYHSGExeqbojWgtUX5bbMs8Q28lEvW7oydXujczpOP9PWOTYlTQeuM7MDwuXbCSaLPYogQ9cZZvYrSV8BbjazK8L9LgYWmNnVKWWeQpjZa3Bw8JD58+e3zf5WWLt2LRMmTOjY+bIEcO8DR86k027bkrak2ZC27w47bceUl/xJ6fPklZ9mTx5ZZXX6fpahrrbVza7Zs2cvyZqkoSgHHDjOSkwK0fL5itDpEMOtgB2AVwGvBOZL2gNQyr6pT5dweqN5ADNmzLA6zRwSp9OzmsyaldX6vmVEa7Xdtl1w4vDRmnnE7T5+7mG89e1/vXlbXgTJ8G23pJ6v7L+RPHvrNktNnLraVle7RhudFvFVwDUWNP9vlbQJmBKu3yW23zQKTBDqDCfL/9zp/N3tPlc7Okkdp1fpdIjhDwhmq0DSPsA44HHgWuAESQOSdgf2Bm7tsG1ODei1YfqO023aGWJ4FfBLYIakVZJOAi4B9gjDDr8NzLGAu4H5wD3Aj4FTPTKlWloVx3YKa6MO0U7hDw+nF2lndMqJZjbVzLYO55K72Mw2mNm7zOwAMzvYzH4W2/9cM9vTzGaY2YJ22dXvtBJb3W6R2/vAXT3axHFK4rlT+pBlS1fmpnzNE+tkitl20M7Y8g37TktdP+6+VUD95v9sN2n13E/X320kjQUWA783s2MkfQ54M7ABeBB4n5k9lVeGi/goo+is8q2KVbvFPC+mvFmyBDzaFhfyIrY4TgV8FLgX2D5cXgjMNbONki4A5gJn5hXguVNGKY2EJ2975KMuOhKyk77kuE1Zn5slT+Qdp2okTQP+CvhGtM7MbjCzjeHizQSRerm4iPcpVQtvp4U8Eu3k5yxaEehmR5PWHe/I7TpfAj4JbMrY/n6gYf+gu1P6mCIulTL+6dHuTx5N19eo32O0XGfVPPnitlz9zMEF975+iqTFsRXzwsGKSDoGWGNmSyTNSh4p6WxgI3Blo7N4S3wUU8YdkjVDT9Fy4uV1kyxbI3+3U4xu38dRwuNmNjP2mhfb9hrgWEnLCcKtD5cUpR2ZAxwDvNMK5EVxER/lFPVtx0n+gMv+oLstAM1cc1G6fW1V4S3t7mJmc8PQ6+kEGVx/ZmbvknQUQUfmsWb2fJGy3J3SJ5QN26sy70gvEY9SyaITYZadOFcn0gQ7pfkKMAAslARBYsAP5B3gLfE+pwpxGC0CHpHVCdrN62yH2LqA1wMzGzKzY8LPe5nZLmZ2UPjKFXBwEe8rkiJUJKqjl6nyuqIUtfEyOy2CVZ7PBXz04O6UPmPB6ovCWeVvqazMfotmiLshOulaic7XrvlKnd7EW+J9Squis2HfacNeR8w+ryLL2keZCJW4SyWtrrrtWomL9miMYXeK4y1xByjXyZXlM46EfNGNcyuzq1Va6bxr1MnZ7X8f8ck02ll+nCx3Urfrop9xEXc2U5W/94jZ59VWyMfdt6r2w+vr7BopMiGHC3pncXeKk0qrP8Q6uVfalTTLSafOD6HRiIu4k0oVP8Q6CHnadfjozfYT5ZvxmZraT9vcKZIuIRg6uiaa7T627Qzgc8CgmT0erpsLnAS8CHzEzH7SLtucfKr80cWFvE4uljq7VboZxtguGkXxFL1Od9WMpJ0t8UuBo5IrJe0CHAmsjK3bj2Do6f7hMf8RJkt3ukCjH0qzLdkjZp+3+dUp4tcSRdLUkXiLNS+PTa8zGq+p27StJW5mN0manrLpiwTpF38YW3cc8G0zWw88JOkB4FCCOTqdLtBocolIyJsRxaf3HGDmyReOWH/S63bmjHD9pAfXb17fagt+weqLhj04ytp8xOzzeMeJu+fOhlQV/SByyU5Qb4W3hgokyWq+8EDEr4vcKZKOBY4ws4+G2btmmtnjkr5CkCMgyuJ1MbDAzK5OKfMU4BSAwcHBQ+bPn982+1th7dq1TJgwodtmpNKMbUVikW38uEJlvTig1PVTJozj8bUbRqwfu97YZ5+dCpWdxW9+80hLx++44wBPPLG+ZTuyaCXWe4edtuPJR56r0JpqqNKuaMRsK8yePXuJmc1spYyd959sH5r/2kL7/sMB17d8viJ0LMRQ0rbA2cAb0zanrEt9uoTpHOcBzJgxw2Z1onnUBENDQ4wm26LdW40lj3h6z4ER60563c5c/D+/Tz/g1odZ/I2PFzp3krLum8i2+L+Bd5y4O9+56iHgoWH7VuXnL1u/cY6fexjXnFfdCNyqqMoub4Hn08k48T2B3YE7wuxc04DbJB0KrAJ2ie07DXi4g7Y5BSn69zfPb75h32nDBDJi7KH5/wpnprhaYIuQxgcbVeF3T3vQPL3nwLDzVx0T38zgJBs/Lveh6dE4o5uOibiZ3Qm8JFpOuFOuBb4l6ULgZcDewK2dss0pRzN+zPj+aaISiVBcINNENG39Zv96KLDt6Dh9cUCbz5vWUq8zRdLrOr1LO0MMrwJmAVMkrQI+bWYXp+1rZndLmg/cQzAl0alm9mK7bHNap4iQp01knDWy7+ipp6J1Ow8L/YtEMkvM06hCYCc9uD7znM/uFnj+Jq5oX19SUcp00Cb37SVR99Gg+bQzOuXEBtunJ5bPBc5tlz1O9eTFM2f92PLWf/3fLgdGxnDniWq7KHLOaHs30gy0GioZHV9EzMvs24jk/e+HaJx247lTnEqIWuattJSCCISgIywZwtioZd0OkU/7J5DWAn96z4ERLpx2iXrVce6NXC3x87XqlsnLBuli3jwu4k5lVDlLUDxhFTQWr3b6pyc9uB5el79PXOgjv3yyg3XcfasK1VE7RmwmH3Lx+ira0m61Jd7MQ96nkGuMi7hTS5I/3qIC0o4RmWVb+dH+M0++kEmx9VHe9aiVnidO8VmXiopYvI7i9VAkyiZ+TLycyLVVlQ+9iJBnzUDlpOMi7tSWZlqkRXOiJFvOeSRb4lHnZhZxl0uaWG5unafYmSaWeal0x923Cq3bGSjv7ohsS9ZFsvyqO0HT+k/6UagljQduIpgYeSvgajP7dLjtNODDBIEe15vZJ7PKcRF3eoKiLdI8AS8astgqz+6mEb7zvHOkuTZGxLzHWsnJFnd8pGyjB1gyuiZpV5qgtyOSpa45bDrMeuBwM1sraWvg55IWANsQpCI50MzWS3pJXiEu4l2k0/MzOp0jLuStPCRmnnwhxMImmxG/tEibtAdNks3hmlQj5HkzQtUpw2WnsCDnydpwcevwZcAHgfPDXFKY2Zq8clzEExQNlav6nC7k7aWq1nYjV0qzJFvARfzY0T7xgUhlOniLX8tAISHPe8A0evjUcWq/NJ7ZOJ6Fj+xbcO/rp0haHFsxL0wbspkwW+sSYC/gIjO7RdI+wOsknQusA84ws19lncVFvAGdEnUX8t4j3pJNE8S81m48fDHPpZEso2jsepot63cdmVwsYmBldvKy4NoCIYfiEUN56QB6bdRrkzzeKAFWOKjxIEmTge9LOoBAl3cAXgW8EpgvaQ/LyFboIu44JcgS5YkrrLCQlx2JmvcweHY3sWmcePx1L+SWMX3aYxy5030j1ketyuUM5h6/fld4drfxgR177jlCfNOuJStbZXzftI7ffsTMnpI0RDCfwirgmlC0b5W0CZgCPJZ2rIt4SapsMftAh+rp5ow9eUIebYfmXDvJcuOt6q0nbGD6tqm/bwCWr8oW6CN3uo+Fj+zL9GmP5e43fdpjQVq6VwfLD64aZMr/bN3Q7iIt7j5plY9A0iDwQijg2wBvAC4g8JMfDgyFrpVxwONZ5biIJygSBVF1LofonO5SyaeVgR+dGvGZJeRQrDOxiJ86LuDTpz0GzwchhvGWdtxvO33aFoF/6/a3AXD1Mwc3PE+SI3e6b8vxOx3Mwmn7snzV4AhXTNao1m6Muq05U4HLQr/4GGC+mV0naRxwiaS7gA3AnCxXCvSRiJcR3jJiEQlvcv/j5x7WkZlgnPK0+y98M0JeRYdpow63hY/sm7pPWgt8YOW4ET70SMCjz2/d/rZUMY+uZdO44ddURKRnnnxh03njew0zWwq8ImX9BuBdRctp68w+7eYV225rvz700Nx9lv5iWcNyDnz13k0fm8WU3bbn8RXPNDxH8nx5+0X2FCkrj6eeeorJkye3VEY7KGpXo/uyadvmWnQvbpM95exOk8bxyNPpHYMbxzd1uqaxgU3DlqeNHc+qF9cxfiDbL75ufWPXR4TWD68HG9jE+IEX2GHc8wxutXbYtsc2BjNEPblh22Hn0PoxTNtmgFV/HPmw3GpdYxvG/nHLNY55Pihj2Pd+aKhxIQkktTzTzsQZO9kh/1FMX//7DV8YXTP71JlWxLobNBL70c6Br947956NeX5900JelqoEPCnMSSHN26ZthNaPYR3ZQh0/JhLRorZHZT/JtiNEPM74gRc2C7kNbMKU3kCMzpsn5i9uM2azkEf38vbbV3LQQa1P0zba6GkRf36XXRo+kc/sUqfh8aclpqa6v7Eb58ypp8L9wee0fePXsmCoed/57TWdOq6MXQfSuEO4lQ7O5F//rKnjqoobT7ou8kL+knxkr535twcC2yJXTdyuIrnPk9eRFj65ftcNDaNcIpavGuQTE3bjC2tXAI2vp1EoZkTd48i7QU+LeBHqFAGS13GZtK/RjODeCdree1tkYoh2UVbwxuxqm9dFojfpwcbnaZRiN77t2d3EwMpxLGeQheH6NDGPmD7tMQbCDlcIHgBlHkxxG/staqUs7ZzZ5xLgGGBNbLb7zwFvJuhxfRB4n5k9FW6bC5wEvAh8xMx+UqU9dUlpmUxNWoQsu33Gk4Cse1tluOGmcapEvLMGCJURuKzO0XiHYFQfedcfj1cvGrseCTnA11cNbo58WR77nEWekBf55+Ct8HTa2RK/FPgKcHls3UJgrpltlHQBMBc4U9J+wAnA/gRzbP5U0j7tmqIta3aRNDHIa+3F928m416nU3yOZvKEHMq7ViIx2xJpsWmYy6OZVmWSpOsjKV55D42kkCdbq0Vm0EmbcCPrH0hWxE08siUZ5bJ+wkh5ieqwbP25gGfTzunZbpI0PbHuhtjizcBbw8/HAd8OE748JOkB4FDgl0XPV1TE8mYXydpedn18qrEk7RyM0u8JtfL+bZWp96f3HBg2TD1wDQwOGxCT1qpsNAw/i7KjQOMs/sbHGRoaYtGN78jdL62RkTVhddZk1XG3SkTeUP488h6IReLpnS00DDEMc94eQ5BR+WXAH4G7CHLc3t3g2OnAdZE7JbHtR8B3zOwKSV8BbjazK8JtFwMLzOzqlONOAU4BGBwcPGT+/PksW7py8/Zgiq+RRPtkba+atWvXMmHChGHn7hSNrjFuW52o0q60Oo+nbM3jhYlj2DQuiIwYGLcRgB03bcsTY54HYP2GoO0zZsMYxmwwNo0TYzaM/B0l46Qjon3Hrg/es4anFynjT3d7aVvu57KlK0fUV9LOuG1RfcV56dgBHn0xeCBE9RgnqseIMRvGxD5vqc+x64199tmphPXpzJ49u/9CDCWdQ+DDHiKY/HANMB7YBzg/FPhPhEHrhZF0NkGy8yujVSm7pT5dwixg8wBmzJhhF5z43cQet4w8KMaC1e8pY2rTDMUiLWbN6mzHavIak+c+86q31TI6ZajCqJmsOm/UGk9rhQP8zfMH8a1tbwdg+RNBa3xg5bjcFmNWKzrPDRLZULSMxXPeUWm9RaTVX7LuorqKE29hR9Epmb7ybYO3pBsmWa/9MvinWRq5U35lZudkbLswTFZeqmkraQ5By/6I2FDSVcAusd2mAQ8XKS/NdRB3rSQ/d4syPvPkzCplaXTNy5au7NvRpHlulTRf8Gb3yYStNos3NBZwKOYOSSPySzcS8E5FbaT51yO3y6QHh4v5MNfIfsXKj3eOQvAgmLhiaxfvgmSPKADM7PoG29eY2eK8feJIOgo4EzjWzJ6PbboWOEHSgKTdgb2BW4uWmyTLx90LvuKkwGzYd1pTPvQoF0uj7XWI2GkXZe53UsDzOt7KdMo149stKuDt7uzLmjZtweqLRjQu0q4z7h4pQrzF3u4wztFEoY5NSTOBs4HdwmNEMDHFgTnHXAXMAqZIWgV8miAaZQBYKAkCP/gHzOxuSfOBewjcLKe2KzKlDjQT7tiuabJgdHeIptV1Xms8GfaX1QFXJpIkT8jLJn3acp6BruYZib4rR8w+b/M/h/g/j4krjDG72uaQxKRLJRljHg0WijqPm+0w7UeKRqdcCfw9cCcwsgcjBTM7MWX1xTn7nwucW9CenidPMOOx5HHSZiOvktEq5mUfmvHIk+h9q11HRkwUbS1WHW0RL2vmyRfy+XeVz0hYFZvnAA1JXmewnC7kcYZlYCQ/fa4znKL/dx4zs2vN7CEzWxG92mpZH7PoxrltE+oi1M3NEnf9NGtX2QdTvCXYrqyDVXHvike7ev5FN84N/ePrR/jpJz24nokrwhZ5CWGePu0xpl92QdWmjkqKtsQ/LekbwCKCGZoBMLNr2mKVA2T/7W+nayVOHQYQpYl2syNVs1rkyYE9eW6TsjR7fNnjZp584WYB7cbAmHiLPLJj7KHxDtgBirTIYUurfCEw/bILWD7nzHaZPSooKuLvA/YlmI05cqcY4CLeJhq5ADol5O0mfo1nXvW2zG15xxcV8iICHqeMkHZ6cEreP4FuzR6fdK3ESQr5QvJzrxBt3+k+Zi16bMsycPb+11Vode9TVMRfbmZ/1lZLnBFEQt6t6cY6TV7oY15/QJnEYvGyWhHwMqKdl0wr7xhI7/jMEvB4DpS4mHZS0PM6jaOQxFWHjy8s5E5jivrEbw7zmzhOpWSFsUXbYHjYZVUPtCiML22C4apa1Wk+4qLHRQ+rsseniX5W67gdNPpXNO6+VUz72bpUH3k0W1DydeRO97nY51BUxF8L3C7pfklLJd0pqdQozX6h6k7BBasv6mryn051cGYNTqraJx89BOKjMht1YjYi3qmX1cFXpqy4nfF1kd1xitp7xOzzNr/azaIb57LoxrloXXaYYNTZ2WhKORg+LdxoQtJ4SbdKukPS3ZL+KVy/o6SFkpaF7zvklVNUxI8iGIDzRoJh+MeE704GVeZL6WRLKo1OR6ukxXXHKdsajwZMbdh3Gk/vOcCqw9OntCkr4GXEuuh+eW6XqiJiOvl9GnffqmGvaN3g9Q8ycYWx+hc7bxbyvMmbR6mQrwcON7OXAwcBR0l6FXAWsMjM9iYIJjkrr5BcEZc0ASAeVpgMMYz2cUaSDI2LxKmMINYp1K/dxBN3NYpAKSrkcf93fCRkq2ll6zBRQStun243DGBL+OHyVYN8ffHrWPjIvrlCPtqwgGi+u63DlxFkdb0sXH8Z8Ja8chp1bP5Q0u3AD4ElZvYcgKQ9gNnA24GvAyOyDTrpNCPK7UxfWzfiUTnD6ioj1DIiLyIjbRh7kTDCuFCPPdQ6JtxFZrOpwm8fCXm3c3VHI2OXr9oyaxDUs/W9YcNWZeLdp0iKpyWZFybw24ykscASYC/gIjO7RdJLzWw1gJmtDnNUZZIr4mZ2hKQ3AX8HvCb0zWwkmAnyemCOmT1S9Ir6lbTkV6N1dGQVNDsL08yTLwS2ZL1La20mh3NP+Z+ta9GqropevZa4kH+9wCxBPcLjjVLRhulFDpI0Gfi+pBFpuxvRMMTQzP4L+K+yBTsBSZGuyzRxvUTev5A0H3Ik5iS2JdPLxltUeTHOVVEm1DDZGt8SZ51ffi8T/3e0nC1i3g+RKWb2lKQhgv7HRyVNDVvhUwlSgGcy6idK7jRlRNpb4dk0Gl2ZpGinX1YL7+ipp7IoEZcf71DVup2HuR2aFfyyMePxfZuJN+8mQR9Hfn7/LCJBLzoxcy8iaRB4IRTwbYA3ABcQZHWdA5wfvv8wrxwX8S5QlXh3atRmXR42aYNziviG4y6UaELfLL9mVt73KmeEKirGyX3yjmmlFd7O9ApZD+P4v6vkXKN9xFTgstAvPgaYb2bXSfolMF/SScBK4G15hbiI15Qyky1D+zIbRrZkTQrdbfLybidn6ImEe/mqwc0tvUkPrsssO36dQ0NDw7Z1wv2SR7LTtc4UnRAl7YG8vjOzKXaFcEa0V6Ss/wNwRNFyCsWJS9pT0kD4eZakj4SOeCdB8kvajP87b9b2LOKx0PFXVcRt6vQAoDhZor1+1w2bX0mSLpRoVp5pPwsEfNx9q0o/mFqtg1b911HIZK+RrOdGA6SaCf/sN4q2xL8HzJS0F0FO8GuBbwFvapdhvUSjH3TW9jJTtsEWIS8jzu0KTWw2k2BZFqy+iCNmn5c7200WRScWaNb+VkM/6+Dj7lYSteR3OarLeM6YSQ+u59nd0gdmOVsoOmJzk5ltBP4P8CUz+xiBPycTSZdIWiPprti6zOGkkuZKeiAc2v+XzVxMt2hWBJrNkd3qD8/Gj2tbK71qIgFPo+gkxctXDW52pRSZG7MIVT24ej2ipBnidZccydlszph+pmhL/AVJJxL0lEbD7bducMylwFeAy2ProuGk50s6K1w+M0yudQKwP/Ay4KeS9umVKdpamaigm+GGVfrT29E5Fm+BT1xhqeK7ZV3n/3bHpyhrhbyMhUWOa5ZupjJO+xca/R7idt1/Y336X+pKmXziHwDONbOHwsmMr8g7wMxukjQ9sfo4gnk3IRhOOkQwcfJxwLfNbD3wkKQHgEOBXxa0r/Y0ytRXVWu8mRZ2VVEuVbtY4qKWnOKsiCsibY7MdswUX1UnZydbn3XKRd9rk5nXDZkV+2sZxjHuamb3Fy48EPHrzOyAcPkpM5sc2/6kme0g6SsEkyZfEa6/GFhgZiOG80s6BTgFYHBw8JD58+cXNacjRImvdthpO5585DmgeHhalUmzbHx2y3THHQd44ol0wcjLPFeWsmF5a9euZcKELal47l49fIzDmA3Bd3Xs+uHf2RcHxKZxjf3l0fFpZeyzz06lbEvjN7/pzuDlvPsZJ+/eVhlCGVGkzjrJ7NmzlzQaQdmIgd2n2dRzTiu074r3ntXy+YpQdLb7NwOfJ/jPurukg4B/NrNjK7Ij7ReY+nQJcw/MA5gxY4bNyppFoEvMmhW0SI+fexjXnHdLqZbFBSdW61rJapW/48Td+c5VD2UeV10rrdz1Dw0NEd3PI2afNyLb4MQVxuD1Dw5bF2UmLJKiNb+l+1BulErctiw+80/dCTlsdD/jZN3bBavfU6VJQLE6c1qnaMfmOQTujacAzOx2YPcmzvdoOIyUxHDSVcAusf2mAQ83UX4tWLD6IvY+cNeu/zVMpgAtSpXhic36/JOukokr0hNQjbtvVaY4R2Frg9c/2NBVUSf3QjtJu6/d/p46rVFUxDea2dOJdc108UfDSWH4cNJrgRMkDYT+9r2BW5so38kgTaSiWOMs33LVseZFmXnyhSNa1s1k9YuuLfL5pz3QogkMqhCyXn0QeC6f3qaoiN8l6W+AsZL2lvTvwC/yDpB0FUHH5AxJq8IhpOcDR0paBhwZLmNmdwPzgXuAHwOn9kpkSi/RTKscaHnwUBmRSBPwSKCzbI9a49F+SUGPhDy+f1pd9IuQ90ta436haHTKacDZBDNRfAv4CfCZvAPM7MSMTanDSc3sXODcgvaMOjrZGqqy87Jq0gb0xOeczGLcfauYxDSyMv09vecAk5jW9dzZdaDsZNNOvSnaEv8rMzvbzF4Zvv4BqKpTs+/pxt/ZPF9yHq20xhtd570rHm2q7Ii0a3p2N21+ZU3LVjW90BpPo9PT8DnVUFTE05ov3qSpgG7+aCLRKyvm7XCtbM4BHiPKb1JWFKPrSWvVz/jMF0uVVYbkSETH6QS57hRJRxPkR9lZ0r/FNm1PMMOPUyFlc6lUQZn8H8kO0MB9UZ60gTGTAF43ct+yYljkemZ85ovc/38/VqrcZqjjtHpF6tNdK71FI5/4w8BiAtfJktj6Z4H2/wpGOWkZD6MfTyeH5BcRm+ToyYD2J29qZjTpuPtWwZ57bk5Hm0bUIo+39KsQruR966aQt/JvwIW8d8h1p5jZHWZ2GbCXmV0We11jZk92yMZRy4LVF21+ZW0vsq7dpAv4cH9z9Equz9svj3iESTMiGLlUopwrRZJeVfXQTN4jd6047aSoT3y6pKsl3SPpt9GrrZb1GXlCHm+dp+2b9yAoSpbQlE3KlBToLMHOm8whTXA7FbM+WoTcHxz9Q9EQw/8EPg18EZhNkBCr7+ZSajd5Qpwm3I32ayapVjuFskwK2KzO1mS8d9nzF/kXUFUir2Rys7i9dfOVO71LURHfxswWSZKZrQDOkfQ/BMLu1IQ0oa9ayNPmQ4wyBUazsJTPWxIw9tD0ofVZJO1sJOpl5+aELQJ8/NzDaDYNSNqDIO++uMD3B5J2IUjVvROwCZhnZl+ObT8D+BwwaGaPZ5VTVMTXSRoDLJP0YeD3wEuaNd7pHK12kGbluo63apMpX+PHJUmb0SWy8+v/dvkIIS7T4cqeew47bzxdbZqAlw2trLKzL6ucKJ+2C3lfsBH4hJndJmkisETSQjO7JxT4IwkmSs6lqIifDmwLfIRgpObhbMmB4owy0qaBSwre03sOjBDysudIClmQNCw7m14yNDHNX5+2LmsUaGRHnahqogmn/pjZamB1+PlZSfcCOxOkH/ki8Em25JfKpFDHppn9yszWmtkqM3ufmR1vZjc3b77TSZrt+MzLghgNEpr2s3XDIkDi2QbjxyfLacaeokPmG0XBtCLgnQr7bCU9QFUPJh+92TJTJC2OvU7J2jGce+EVwC2SjgV+b2Z3FDlJ0Xzi+wB/D+wWP8bMDi9yvNP75P3F7+SMNHFxO3rqqZttimbXSWuJR/8Y4n7wVmOoI+oYS13GFZNXDxv2nTZs6jQHtEFl/nU+XmRSCEkTCCajP53AxXI28MaiJynqTvku8DXg64BnF+xR4v7x5I+yUaurjj7apL8/EvisIfzNtsAXrL6IoaEhFqx+T+4ArV6k0TyrdbzvowlJWxMI+JVmdo2kPyOYq+EOSRDMrXCbpEPNLHXqqDL5xL9qZrea2ZLoVcVFOJ0lmrCiDNEPuVH+8TJUGY+dFNHF3/j4ZndP/BVRRpiKhna2gyjXeSdoVCfuWqkeBSp9MXCvmV0IYGZ3mtlLzGy6mU0nmDDn4CwBh+Ii/iNJH5I0VdKO0avVi3B6lzJC3o2W6qIb57bsGy46krbdAlcnIXcxr5TXAO8GDpd0e/h6U9lCior4HAKf+C8IcqgsIcip4vQxjYQ8bbRpp2hVbBrZW8Uo2TJ0Iw96Vqe2C3k1mNnPzUxmdqCZHRS+/iuxz/S8GHEoHp2ye8prj2aNl/QxSXdLukvSVZLGh637hZKWhe87NFu+Uy3NtGi7OflCJDJZLcsovLHXfNntdq/kzZyUxIW8PuSKuKTDw/fj017NnFDSzgTx5jPN7ABgLHACcBawyMz2BhaFy04HKPKDzJpEogr/eNXkiXM7pmSDzopau4S8bCdm5F5xQe8ujVrirw/f35zyOqaF824FbCNpK4JBRA8DxwGXhdsvA97SQvlOCcq0Susi5I2EI3k98QmRi3RWNiNMnRbyIlPWNTuvqtM7yKyZSetbPKn0UYL5NP8I3GBm75T0lJlNju3zpJmNcKmEAfOnAAwODh4yf/78DlldjrVr1zJhwoRum5FKnm3LlmaP8rXxQXzsiwPpg2jGrt/yXdpnn50alp2MkilTZ1FZeZE2v/nN8A79NJsa2RidI8u2rP07QXTuHXbajicfeW7E+ZO2Rfcvj1bmX23lfnaC2bNnLykSt53H+J13sd0+8PFC+/7mHz/e8vmK0Ghmn1xro7CYMoS+7uMIYiGfAr4r6V1FjzezecA8gBkzZtisZrMStZmhoSF60bZZs7JblPG/22mt73grfdGNJ4zYfsGJw8tNDrEvU2dRWVnD9JPD1tOG+eeVm7Qzy7as/TtBdK+On3sY15x3S3jui0Zsj9PIZdJKq72V++k0TyN3ysTwNRP4IMG4/p2BDwD7NXnONwAPmdljZvYCcA3wauBRSVMBwvc1TZbvtIFGP/4iozbb0ZGY9sBJrisqTMmRmL3Q+Rm3r4iteXXRqtvFfePdIbclbmb/BCDpBoKA82fD5XMIRnE2w0rgVZK2JXCnHEEQrvgcQSjj+eF7w8QvTntIy3yYHHYfZQiMi3fUOs8T9E4JeZJ2jTxMq6tOj+JslDgsSTuzJPb6CNZepGic+K5A3Fm2AZjezAnN7BbgauA24M7QhnkE4n2kpGUEKRjPb6Z8pxqK/BCzxLpOEStR63LcfatKdYaWaVWmtdh7rVXqnZ+9S1ER/yZwq6RzJH0auIUtkSSlMbNPm9m+ZnaAmb3bzNab2R/M7Agz2zt8f6LZ8h2n2YEqSSFvRYzrIORZD+Os7JJO79FQxMPx/ZcTTMn2JEFn5PvMzBMe9xmj5cdeVlyL7F8HwU6jrnY51dEwi6GZmaQfmNkhBC4QxxnVlJ0NKW3fOviFOy3gdbjmfqRoKtqbJb3SzH7VVmucWpHm583qFKuTH7wKkkK+bOnKzDk26yheWTnPvWU++ijqE59NIOQPSloq6U5JS9tpmFM/ssQqTcBHw/RivdZZmTYMvpupdJ3OUFTEjwb2IJhbMxpy/+Z2GeX0Np2c6adZigpyL8SKQzmXTi9cj1OcolkMVwCT2ZI3ZXK4znF6ljLRJ0Esdj3FLy9VQhbteDjV/Z/KaKWQiIe5Tq4EXhK+rpB0WjsNc3qXyL3iP+r2keU66fY/B7/nnadox+ZJwGFm9hyApAuAXwL/3i7DnHqyYPVFI/zd0ejNNHwEX3XkCWS3fd9+j7tHUREXwydIfjFc5zhOh0kK5tDQUEtlVdF69od19ygq4v8J3CLp+wTifRzBBJ9OH5IWZpjWGt+w77S2DRDqx7A5F8nRhaRLCIJE1oQT5CDpIOBrwHhgI/AhM7s1r5yiHZsXEozYfAL4A8GIzS81a7zT2yxYfVHh2WXalWgpTqvi1i8PgSy66UeP+/b78D5cChyVWPdZ4J/M7CDgH8PlXIq2xCFwoVj42lTiOMepPWXcAVnhfHUduVlXkvXVb3VlZjdJmp5cDWwffp5EMOtZLmWjU6bg0SnOKKXZHCmR+FQ1zVuvEr/Wo6eemhr6mNXqHqUCPkXS4tjrlALHnA58TtLvgM8DDf/yenSK0zTtzEtdlqo66PIoUn4kRklBG6UiNYK03OqN6KW6GbsBJq4oPKXl401Mz/ZB4GNm9j1Jbyfoe3xD3gEeneK0nU6JWBVCHh2ftDetVZl3Tb0kTN3A6yeTOcBHw8/fBb7R6IBmolMgmIneo1Oc3IRY8eH3WeJYV5IPg+PnHjZsuVeuo254vTXkYeD1wBBBmpNljQ4oJOJmdqGkIeC1BC3w95nZr5u1UtJkgifMAQSO/PcD9wPfIZgxaDnwdjN7stlzOP1JO90qLkCt4fU3HElXAbMIfOergE8Dfwt8WdJWwDqgoR+9kIhLehVwt5ndFi5PlHRYONVaM3wZ+LGZvVXSOGBb4FPAIjM7X9JZwFnAmU2W73SBvJGbvcyC1Re1NKCmV0jz51ddtrMFMzsxY9MhZcopmsXwq8Da2PJz4brSSNoe+AtCd4yZbTCzpwgGEEVTvl1G4LJxak7clVIHAa9SgLqdh6RbNHPN3j/QPWTWuKdV0u1h8Hl83VIzO7D0CYMRSfOAe4CXA0sIHPm/N7PJsf2eNLMdUo4/hfAvxuDg4CHz588va0JHWLt2LRMmTOi2GalUZduypSux8eM2L784MLyve+z64LuldVvm2N77wF3balczGf3SSNrZD/czSZm6jOorfszUPXasVZ3Nnj17SRPRIsPYbnAX+9PjPlZo3yUXf6Ll8xWhaMfmbyV9hC2t7w8Bv23hnAcDp5nZLZK+TOA6KYSZzSN4CDBjxgyblTXdSpcZGhpiNNt29NRTR3RoxlviUafmyGH3t2S2zKqw64ITq2qJD7dztN/PNIrWZbyeomMiF1Rd62w0UdSd8gHg1cDvgVXAYRRwuGewClgV86dfTSDqj0qaChC+r2myfKfN5LksJj24vuGkEP00AGa0k5Y90d0nnaVodMoa4IQqTmhmj0j6naQZZnY/cASBa+UeghjJ88P3H1ZxPqezRKGFjRJf1XEATD8m1XJ6nzK5UwCQdJuZHdzieU8DrgwjU35LkFxrDDBf0knASuBtLZ7D6QK9MDVbGt3Ox91N/IHV2+SKuKT/IkiFuDy+utWTmtntQJrD/4hWy3a6Q3JwTxZ1FMc62tQNOpG6wKmeRi3xS4EbJF0GfNbMXgCub7tVTk9Rh9BCp3niD7F2xoo77SG3Y9PM5gOvIEiNuFjSGcATkj4u6eOdMNCpP73iQvFOt+pxse8+RaJTXiAY3DMATEy8HAeol5A3EurkhMIuRCMp87Dz+usujXziRwEXAtcCB5vZ8x2xyqkVyUiSZn2nncpk2My56hgt003K3l+vv+7RqCV+NvA2MzvLBbw/iX7MyfzYSZ7ec2CYbzw5EMh/4P1Bn06z1lVyW+Jm9rpOGeLUn2Z/nHUW8Drb1k3K/ttq9h+Q0zpFR2w6zgjiA3ryRmp6y2z04p3F3af0YB+nf4iLb9ZoxvikEGm5U5zRi4t3PfCWuNMQ/7H2H8kHeD+PaK073hJ3MvEfav/RJ7PQjypcxJ3SeHz16CLrPnpCsN7A3SlOR3AR6B3yOiu9ZV4dki6RtEbSXbF1n5N0n6Slkr4fzkeci4u4UzkepdKbNJpizcMIK+dS4KjEuoXAAeGsab8B5jYqxEXcqYxGE0L4j79+RA9WF+jOY2Y3AU8k1t1gZhvDxZuBaSMOTOAi7pSmyKi8tEkhvCXu9BlTJC2OvcrOhvZ+YEGjnbxj02kLG/ad1nB2H6f7eMu7HGPXWZkxEI83O1GypLOBjcCVjfbtWktc0lhJv5Z0Xbi8o6SFkpaF7yNmund6GxeMeuL5TuqFpDnAMcA7zcwa7d9Nd8pHgXtjy2cBi8xsb2BRuOw4TodwMe8+YebYM4FjiyYd7IqIS5oG/BXwjdjq44DLws+XAW/psFlOk2S5TZKZDCORcLGoF/4PqTtIugr4JTBD0qpwfuGvEMzVsFDS7ZK+1rCcAq31ypF0NXAegbFnmNkxkp4ys8mxfZ40sxEulbBz4BSAwcHBQ+bPn98hq8uxdu1aJkyY0G0zUmnVtmVLV45YZ+PHpe6rdRsyy9n7wF0rtaud9INt0X1N3pdmqVudzZ49e0mzPuqI7SdOs1fOLNYA+dnQp1o+XxE63rEp6RhgjZktkTSr7PFmNg+YBzBjxgybNat0ER1haGiI0WrbBSeO/BInW92NOjXTWn+juc7aSRW2Dc+V8p4WLQqoc52NJroRnfIa4FhJbwLGA9tLugJ4VNJUM1staSqwpgu2OU3ikSi9ibu1ep+O+8TNbK6ZTTOz6cAJwM/M7F0EU8DNCXebA/yw07Y5ncF9sPXABXx0UKfBPucDR0paBhwZLjs1pBURdgGvBy7go4euDvYxsyFgKPz8B+CIbtrjFMMFoLcpkrXQ6R3q1BJ3+gB/AHSPZGinJ7UaHbiIOx3Hhbzz5E324ALe27iIO13BhbxzNJpqzeltXMQdp09w8R6deBZDxxnluHiPbrwl7nQFFxbHqQYXcacU7st2nHrhIu44jtPDuIg7hamyFZ6WknbZ0pXe0neckriIO4Vop7h6fnHHaR4XcachLrCOU19cxJ1cXMAdp954nLiTSTsFPGtSiKomJHCcfsFb4o7jOF1C0mRJV0u6T9K9kv68bBneEnccx+keXwZ+bGZvlTQO2LZsAS7ijuM4XUDS9sBfAO8FMLMNQPbM4hl0Y6LkXYDLgZ2ATcA8M/uypB2B7wDTgeXA283syU7b5wTk+cOzhswX8aH7cHvH2cwewGPAf0p6ObAE+KiZPVemkG60xDcCnzCz2yRNBJZIWkjwNFpkZudLOgs4CzizC/b1JUU7MV2EnX5G6zaUmRR8iqTFseV5ZjYvtrwVcDBwmpndIunLBLr3f8vY1HERN7PVwOrw87OS7gV2Bo4DZoW7XUYwbZuLeI1wAXecUjxuZjNztq8CVpnZLeHy1QQiXgqZWTPGVYKk6cBNwAHASjObHNv2pJntkHLMKcApAIODg4fMnz+/M8aWZO3atUyYMKHbZqSStG3Z0pWFjtv7wF0b7tOorLwyeqnO6kRdbaubXbNnz17SQFQbMmnrl9irp7yt0L4/fuQ/Gp5P0v8AJ5vZ/ZLOAbYzs78vY1PXOjYlTQC+B5xuZs9IKnRc+HdkHsCMGTNs1qxZbbOxFYaGhugF26r2Y8+a1cifnh0H3it1Vjfqaltd7aoZpwFXhpEpvwXeV7aArsSJS9qaQMCvNLNrwtWPSpoabp8KrOmGbf1EN0Zj+ghQx9mCmd1uZjPN7EAze0szwRwdF3EFTe6LgXvN7MLYpmuBOeHnOcAPO21bPxPNvdgJv7cLueNURzfcKa8B3g3cKen2cN2ngPOB+ZJOAlYCxRxPTlMkJ8+NU6WQJ8tyAXecaulGdMrPgSwH+BGdtMVpb8RJWtke4eI41eIjNvuUTou34zjtwRNg9SE+g47jjB5cxPsMF2/HGV24iPcZ7upwnNGFi3gfEo2a9Fa54/Q+LuKO4zg9jIt4HxLlN3HXiuP0Pi7ifUbcheLuFMfpfTxOvE/xVrjjjA68Jd5HRC3vIillHcfpDVzE+5Ci+cMdx6k/LuJ9RNyF4v5wxxkduIj3Ke4Td5zRgXdsjlLyWtruE3ec0YOLeB8Rtb6Hhoa6a4jjOJXhIj4KiVrh7jJxnNFP7Xziko6SdL+kBySd1W17eg3vsHSc3qAqrauViEsaC1wEHA3sB5woab/uWtVbeOvbcepPlVpXN3fKocADZvZbAEnfBo4D7umqVT2GC7nj1J7KtK5uIr4z8LvY8irgsPgOkk4BTgkX10u6q0O2lWUK8Hi3jcigrrbV1S5w25qhbnbt1moBz2x87Cc/fuQ/phTcfbykxbHleWY2L/zcUOuKUjcRT5tA2YYtBJUwD0DSYjOb2QnDyuK2laeudoHb1gx1tasVzOyoiopqqHVFqZVPnOBptEtseRrwcJdscRzHaReVaV3dRPxXwN6Sdpc0DjgBuLbLNjmO41RNZVpXK3eKmW2U9GHgJ8BY4BIzuzvnkHk527qN21aeutoFblsz1NWurtOE1mUis6bcMI7jOE4NqJs7xXEcxymBi7jjOE4P07MiXpfh+ZJ2kXSjpHsl3S3po+H6cyT9XtLt4etNXbJvuaQ7QxsWh+t2lLRQ0rLwfYcu2DUjVje3S3pG0undqDdJl0haEx9zkFdHkuaG37v7Jf1lF2z7nKT7JC2V9H1Jk8P10yX9MVZ3X+uCbZn3r5P11leYWc+9CDoCHgT2AMYBdwD7dcmWqcDB4eeJwG8IhtGeA5xRg7paDkxJrPsscFb4+Szgghrcz0cIBmN0vN6AvwAOBu5qVEfhvb0DGAB2D7+HYzts2xuBrcLPF8Rsmx7fr0v1lnr/Ol1v/fTq1Zb45iGrZrYBiIasdhwzW21mt4WfnwXuJRiNVWeOAy4LP18GvKV7pgBwBPCgma3oxsnN7CbgicTqrDo6Dvi2ma03s4eABwi+jx2zzcxuMLON4eLNBDHGHSej3rLoaL31E70q4mlDVrsunJKmA68AbglXfTj8y3tJN1wWIQbcIGlJmLIA4KVmthqChxDwki7ZFnECcFVsuQ71llVHdfvuvR9YEFveXdKvJf23pNd1yaa0+1e3ehs19KqIVzZktSokTQC+B5xuZs8AXwX2BA4CVgNf6JJprzGzgwmypZ0q6S+6ZEcq4UCHY4HvhqvqUm9Z1Oa7J+lsYCNwZbhqNbCrmb0C+DjwLUnbd9isrPtXm3obbfSqiNdqeL6krQkE/EozuwbAzB41sxfNbBPwdbr019HMHg7f1wDfD+14VNLU0PapwJpu2BZyNHCbmT0K9ak3suuoFt89SXOAY4B3Wuh0Dl0Vfwg/LyHwO+/TSbty7l8t6m000qsiXpvh+ZIEXAzca2YXxtZPje32f4COZ1uUtJ2kidFngg6xuwjqak642xzgh522LcaJxFwpdai3kKw6uhY4QdKApN2BvYFbO2mYpKOAM4Fjzez52PpBBXmqkbRHaNtvO2xb1v3rer2NWrrds9rsC3gTQSTIg8DZXbTjtQR/C5cCt4evNwHfBO4M118LTO2CbXsQRATcAdwd1RPwJ8AiYFn4vmOX6m5b4A/ApNi6jtcbwUNkNfACQYvxpLw6As4Ov3f3A0d3wbYHCPzL0ffta+G+fx3e5zuA24A3d8G2zPvXyXrrp5cPu3ccx+lhetWd4jiO4+Ai7jiO09O4iDuO4/QwLuKO4zg9jIu44zhOD+Mi7jQkzNT4kKQdw+UdwuWWZw9v0p614fvLJF3dQjmnS9q2IpveIukfSx7z0y6mFXBGCR5i6BRC0ieBvczsFEn/D1huZud14Lxb2ZZkT9G6tWY2oYKylwMzzezxCsr6BcHgm8JlhaMup5nZua2e3+lfvCXuFOWLwKsknU4wwCk1p4mk94TJj+6Q9M1w3W6SFoXrF0natcH6SyVdKOlG4IJwZO4vJf1K0mdi55oe5bKW9F5J10j6sYIc4J+N7fdVSYsV5Hv/p3DdR4CXATeG50HSG8Pz3Cbpu2E+HCSdL+me0M7Pp1zzPsD6SMBD+7+qIM/8byW9PkwGda+kS2OHXkswYtVxmqfbo4381Tsv4C8JRqcembF9f4LReFPC5R3D9x8Bc8LP7wd+0GD9pcB1hPmmCcTuPeHnU4G14efphLmsgfcSDDGfBIwHVgC7JOwYCwwBB4bLy2O2TgFuArYLl88E/hHYMbym6F/r5JTrfh/whdjypQTpkUWQgvUZ4M8IGk1LgINi+y4D/qTb99ZfvfvylrhThqMJhlkfkLH9cOBqC1ukZhblmv5z4Fvh528StOTz1gN818xeDD+/hi35Vb6ZY98iM3vazNYB9xBMMgHwdkm3Ab8meNDsl3Lsq8L1/yvpdoJ8KbsRCPA64BuSjgeeTzl2KvBYYt2PzMwIhqA/amZ3WpAU6m6Ch0/EGoJ/BI7TFFt12wCnN5B0EHAkgdj9XNK3Cb4/Pwp3+RpBy7NIJ0vWPvH1zxU8Js762OcXga3CZEtnAK80sydDd8b4lGMFLDSzEe4NSYcSTFxxAvBhgodVnD8S/ANIs2VTwq5NDP/djQ+Pd5ym8Ja405AwU+NXCXKlrwQ+B3zezH5nZgeFr68RJIp6u6Q/CY/bMSziFwQCCPBO4OcN1if538R+Zdie4IHwtKSXEvybiHiWYEo9CGbIeY2kvULbt5W0T+gXn2Rm/wWcTpAnO8m9wF4l7YrqdScCt47jNIWLuFOEvwVWmtnCcPk/gH0lvT6+k5ndDZwL/LekO4AoNe9HgPdJWgq8G/hog/VJPkowocWvGNnizcXM7iBwo9wNXELwQIiYByyQdKOZPUbgV78qtOdmYF8Ckb8uXPffwMdSTnMT8IpQlMtwCHCzJaJvHKcMHmLoOBUg6csEfvCfljzmWjNb1D7LnNGOt8Qdpxr+lSA/ehnucgF3WsVb4o7jOD2Mt8Qdx3F6GBdxx3GcHsZF3HEcp4dxEXccx+lhXMQdx3F6mP8PZaqwibqI1uoAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "## Extract coords of location, depth, and temperature at a particular location\n",
    "## For illustrative purpose, take point x == 140\n",
    "pt = 125\n",
    "plt.contourf( water_depth)\n",
    "plt.grid()\n",
    "plt.xlabel('X-coordinates (m)')\n",
    "plt.ylabel('Y-coordinates (m)')\n",
    "plt.colorbar()\n",
    "plt.hlines(pt, xmin = 0, xmax=(water_depth.shape[1]), colors= 'r')\n",
    "plt.title('Vertical section through red line')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "<ipython-input-3-268344256d76>:12: DeprecationWarning: `np.bool` is a deprecated alias for the builtin `bool`. To silence this warning, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.\n",
      "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n",
      "  var0 = ncdataset.variables[variable][:]\n"
     ]
    }
   ],
   "source": [
    "loc = x[pt]\n",
    "d_pt = water_depth[pt, :]\n",
    "tracer = getvariablevalues(ncdataset, 'Tracer') ## We extract Tracer to illustrate\n",
    "sigma_layers = getvariablevalues(ncdataset, 'Depth')  ## the number of layers and thickness of each\n",
    "\n",
    "## We extract tracer values corresponding the the vertical line at the pt highlighted in figure above\n",
    "trace_sect = tracer[:,pt,:] # note that tracer is 3D of size n_sigma_layers, X, Y\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Extract the x and y coordinates to plot vertical profile\n",
    "The `tracer_sect` variable generated above is a 2D array of dimensions [nlayers, n_y_gridcells] <br>\n",
    "We also need x and y coordinates to plot <br>\n",
    "We extract y coord for each point in the horizontal and vertical <br>\n",
    "i.e. create an array of size [nlayers, n_y_gridcells] where we specify position in the \n",
    "vertical water column (based on water depth and layer number) and the point in the horizontal (based on the x-coordinate)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "d_pt_arr = np.zeros([len(sigma_layers), len(d_pt)]) # create empty array\n",
    "for i in range(0, len(d_pt)):\n",
    "    # for each grid cell, create the position at the different layers\n",
    "    # note that water depth is positive so we flip\n",
    "    # Also note that sigma layers is from 0 - 100 so we convert to 0 - 1\n",
    "    d_pt_arr[:, i] = -d_pt[i]*(sigma_layers/100) \n",
    "\n",
    "x_plot = np.tile(np.arange(0,len(x)), (len(sigma_layers),1))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Now plot the vertical profile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAEKCAYAAAAcgp5RAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAk/ElEQVR4nO3dfZgddXn/8feHPBISQMxigCRNwLApsTQCUpWiGx4UKQZ8bEQqPrQpLVK15apJ0wrqlYpI689W+8Mo/KRVxEh4KqiQIEdtFTCBEBKTSBICCSRABI0rZPN0//6Y2XCyOWf37Hmas7Of13Wda2e+M2fmzpzNfWa/8517FBGYmVl+HZR1AGZm1lhO9GZmOedEb2aWc070ZmY550RvZpZzTvRmZjnXsole0jmS1kpaJ2lO1vGYmTWCpOslPStpZVHbEZIWS3os/fmKHu+ZKKlT0uWV7KMlE72kIcBXgLcBJwDvk3RCtlGZmTXEN4BzerTNAe6NiCnAvel8sS8C3690By2Z6IFTgXURsSEidgI3AednHJOZWd1FxI+B53s0nw/ckE7fAFzQvUDSBcAGYFWl+xhaU4SNcwywqWh+M/BHxStImg3MBhg5cuTJEydO7HOjES/1uY50cE3vL//eoUi7K9pPVvbu3ctBB7Xqd3+iPzHW8nnVoudn3WqaGd/uKs8lD4oh7OLAu/aHaO++6aHs3bePPXEQsUcogL0BByn5mRo5cti+6V/+8pfbIqKtqsBSp3eMjBee39v3isCqR3etAnYUNS2IiAV9vO1VEbEFICK2SDoSQNIhwCeBs4GKum2gdRO9SrTt96mnB2oBQHt7e6xdu7bPja7ddHSf67RPKL+dSt5fzsY1lzNp6lUV7ScrhUKBjo6OrMPoVX9irOXzqkXPz7rVNDO+m7efVNX7jtt0AdcOW3NA+9njXm5796EP7dvH4q1T2fLTYxj/wySf/ua4ERy2vmvfuvfeN3fftKQnqgqqyAvP72XRXWMrWnfqxC07IuKUWveZ+jTwxYjolEqlydJaNdFvBiYUzY8Hns4oFjOzZntG0lHp2fxRwLNp+x8B75Z0NXA4sFfSjoj4cm8ba9VE/3NgiqTJwFPALODCbEMyM2uaO4CLgavSn7cDRMTp3StIuhLo7CvJQ4tejI2I3cBHgbuB1cDCiKj4woOZDS4bN9fU5Z4pSd8Gfga0S9os6SMkCf5sSY+R9MfX1NfWqmf0RMT3gO9lHYeZWSNFxPvKLDqzj/ddWek+WvKM3szM6seJ3swaotoRN1Z/TvRmZjnnRG9mlnNO9GZmOedEb2aWc070ZmY550RvZpZzTvQZaJ/gsj1m1jxO9E3mJN98PuY22LVsCYQ8csJpLh/v7PhmqdYyqBJ9rf/x+/P+nrXQRww7saZ9W23aJzydWX16s6wNqkTfTD2/FLasL2QTiNkAtH3XSBjW93pWGffRm5nlnBO9DRrus88XXweonBO9mQ1oi7dOBWDMEwc+TNwSLZfoJV0p6SlJy9PXuVnHZGatyWf1lWm5RJ/6YkRMT19+ypTVjbtvBp9775ubdQiZa9VEb2Y2KEj6mKSVklZJ+nhR+2WS1qbtV9eyj1YdXvlRSR8AlgJ/FxEvZB2QmVm9SXoN8BfAqcBO4AeS7gLGA+cDJ0ZEl6Qja9pPRPMvYEhaAowrsWgecD+wDQjgs8BREfHhEtuYDcwGaGtrO3nhwoWNC7gOOjs7GT16dNZh9Gowxdi1a0Udoimz7R3jGDFya8O2X6tmxPfCnlE1vX9o1xE8f9CLB7QfOmzHAW3bd42ka+dQhnaKIV1JPtsz4uXp44/fP9XMmDFjWUScUkt8rzlxeCy6a2xF606duKXs/iS9B3hrRPx5Ov9PQBdwCrAgIpbUEme3TM7oI+KsStaT9DXgzjLbWAAsAGhvb4+Ojo66xdcIhUIBx1i7esW4dtOFtQdTxsY1c5g09aqGbb9WzYhvaY0XSds2zOLGUcsPaD973JoD2hZvncrG59sY+5NhHLa+C4DfHDdi3/S9982qKZYGWwnMl/RK4CXgXJKejOOB0yXNB3YAl0fEz6vdSct13Ug6KiK2pLPvIDkQZnXlkghWixf2jOrHiJ+7xkpaWtSwID1RJSJWS/o8sBjoBB4BdpPk5lcArwdeByyUdGxU2QXTcokeuFrSdJKum43AX2YajZlZbbb11lUUEdcB1wFI+mdgM/D7wC1pYn9Q0l5gLPBcNQG0XKKPiD/LOgYzy4eBMLRS0pER8aykicA7gTcAe4EzgIKk44HhJNcuq9Jyid6sWdx9Yy1iUdpHvwu4NCJekHQ9cL2klSSjcS6uttsGnOjNzDIVEaeXaNsJXFSvffiGKTOznHOit0HNJRFa28bNbVmHkAtO9GZmOedEb2Z1Ve+Kkj6rr50TvZnl0kAYWtksTvRmZjnnRG+Dni/IWt450ZuZ5ZwTvZnlTnflSks40ZuZ5ZwTvZlZzjnRm+ELspZvTvRmZjnnRG9mlnNO9GZmOZdJopf0HkmrJO2VdEqPZXMlrZO0VtJbs4jPBif301teZfXgkZUkj8z6anGjpBOAWcA04GhgiaTjI2JP80M0M8uHTM7oI2J1RKwtseh84KaI6IqIx4F1wKnNjc7MBjoXNNtfqz1K8Bjg/qL5zWnbASTNBmYDtLW1USgUGh5cLTo7Ox1jHTQ6xq5dc2rfxo5xbFxT+3YapZHxvbBnFMfVYTtD947iwhen0zU6SVEjXkzSwKGbph6wbtuukewaPpyDTg+GnJo8VrXVf4+brWGJXtISYFyJRfMi4vZybyvRVvKBuBGxAFgA0N7eHh0dHdWE2TSFQgHHWLvGx9hxwAPDi/vuK3mY+MY1c5g09aq6R1YvjYxvaZ1q0bdtmMWNo5az8fmkFv2k8c8BcPa4NQesu3jrVLYsP4YxT8S+0gf33jerLnE0g6RPAH9OkuseBT4ETAWuBUYCu4G/jogHq91HwxJ9RJxVxds2AxOK5scDvkJmZrkk6Rjgb4ATIuIlSQtJrlNeCHw6Ir4v6VzgaqCj2v202vDKO4BZkkZImgxMAar+FjMzGwCGAgdLGgqMIjm5DeDQdPlh1HjCm0kfvaR3AP8OtAF3SVoeEW+NiFXpN9ovSP5cudQjbqzZ2ic8XbaLprdlPdeDyrp6bODZvnski7ceeL2gtLvGSlpa1LAg7XomIp6SdA3wJPAScE9E3CNpE3B3uuwg4I21xJtJoo+IW4FbyyybD8xvbkRm++tO6B5bb3WwLSJOKbVA0itIRhtOBn4NfFfSRSSjDT8REYskvRe4DqimOxxova4bs5bhJG9NcBbweEQ8FxG7gFtIzt4vTqcBvkuNw8yd6M3MsvMk8HpJoyQJOBNYTdIn/+Z0nTOAx2rZSauNozfLlUr79G1wiogHJN0MPERyXfJhkmHjDwNfSi/Q7iC9Z6haTvRm/dRX8h4x7MQmRmMDXURcAVzRo/l/gJPrtQ933ZiZ5ZwTvVkVfKG2dQ1fsznrEFqOE72ZWc450ZtVqdRZfaVtVj8bN7dlHULLc6I3s1zoLmj2/S1fyTiS1uNEb2aWc070ZjUo7pZxF03z9aw34weOlOZEb1YHfSV5fwk0RuHMa7IOYUDwDVNmNXISz9bGiz+ZTPxTtnG0Mp/Rm5nlnBO9mdXFzXV6jKDVnxO9mbWUyh/oYZXKJNFLeo+kVZL2SjqlqH2SpJckLU9f12YRn5m1vnnT7sw6hAEjq4uxK4F3Al8tsWx9RExvbjhmZvmV1aMEVwMkdfbNzKyRWnF45WRJDwPbgX+MiJ+UWknSbNJi/G1tbRQKheZFWIXOzk7HWAcDO8Yb6dq1otnhHKBrxzg2rplT9+0et2dUXbbTtmskR+wdxYUvTqdrdJKiRrx4DACHbkr678eNfHXL/x60koYleklLgHElFs2LiNvLvG0LMDEifiXpZOA2SdMiYnvPFdOnqC8AaG9vj46OjjpF3hiFQgHHWLuBHuPaTRc2N5gSNq6Zw6SpV9V9u0vrNOpm8dapXPjidG4ctZyNzycFyyaNf27fct8k1X8NS/QR0e8nlkdEF9CVTi+TtB44Hlha5/DMrMW5KmX9tNTwSkltkoak08cCU4AN2UZl1hp8B24+SfpEOgpxpaRvSxop6QhJiyU9lv58RS37yGp45TskbQbeANwl6e500ZuAFZIeAW4GLomI57OI0cys0SQdA/wNcEpEvAYYAswC5gD3RsQU4N50vmpZjbq5Fbi1RPsiYFHzIzIzy8xQ4GBJu4BRwNPAXKAjXX4DUAA+WcsOzCwD7ROeZu2moxv+njxplQuxO3cO7c81hLGSiq8zLkgHkxART0m6BngSeAm4JyLukfSqiNiSrrNF0pG1xOtEb9ZE7mcflLZFxCmlFqR97+cDk4FfA9+VdFG9A2ipi7Fmg40T/6B3FvB4RDwXEbuAW4A3As9IOgog/flsLTupONFLOqR7RIyZWTPleKjlk8DrJY1SUirgTGA1cAdwcbrOxUC5e48qUrbrRtJBJFd/3w+8jmR8+whJzwHfI+lneqyWnZuZDWYR8YCkm4GHgN3AwyQ3go4GFkr6CMmXwXtq2U9vffT3AUtIrv6ujIi9AJKOAGYAV0m6NSK+WUsAZta3Vu/iaUYt+la5EFtvEXEFcEWP5i6Ss/u66C3Rn5X2GfUM6nmSIZCLJA2rVyBmg1V/R9IMxpE3+x4XaFUpm+iLk3x6ZXhC8foR8VCpLwIzs3pxgq+PPodXSvos8EFgPRBpcwBnNC4sMzOrl0rG0b8XOC4idjY6GLPBajB2x1jzVDK8ciVweIPjMDPbj7tt6qeSM/rPAQ9LWklaQhggImY2LCozM6ubShL9DcDngUeBvY0Nx8xaVfcQT3cxDTyVJPptEfFvDY/EzMwaopI++mWSPifpDZJO6n41PDIzK6vVb6Cq1fFjxmcdQq5Uckb/2vTn64vaPLzSzGyA6DPRR8SMeu9U0heAtwM7Scbnfygifp0umwt8BNgD/E1E3F1uO2bWHMV/QfQcCtqM8gdWm7JdN5IuSgublVt+nKQ/rnK/i4HXRMSJwC9J6ukg6QSSQmrTgHOA/3DFTBvs8t5NY43XWx/9K0mGVV4v6VJJ75X0AUmfkfQj4GrgmWp2GhH3RMTudPZ+oLtD7nzgpojoiojHgXXAqdXsw2ygyUtCnzftzqxDsB4UEeUXJmfTZwCnAUeRPOpqNfD9iHiyLgFI/w18JyK+KenLwP3dFTElXZfu6+YS75sNzAZoa2s7eeHChfUIp2E6OzsZPXp01mH0yjHWRy0xdu1acUDbiGEnVrxuRfvYMY4RI7f26z09Yyje9wt7Ru23bNzIV7N1x7qqYtu+ayQAR+vwhn3OM2bMWFbuiU+VGjF5fBx15WUVrfvEB+fUvL9a9dpHHxF7SLpZFvd3w5KWAONKLJoXEben68wjqcH8re63lQqjTGwLSOo2097eHh0dHf0NsakKhQKOsXZ5j3HtpgsPaCt3pl9q3UpsXDOHSVOv6td7DoyhY18//dIeffSzpt3J/FXVlRRevHUqhTOvGRCf80DSsGfGRsRZvS2XdDFwHnBmvPxnxWaSKpndxpM8Ed3MMtKfLiV327SmTJ4ZK+kc4JPAzIh4sWjRHcAsSSMkTQamAA9mEaOZWV407Iy+D18GRgCLk8ckcn9EXBIRqyQtBH5B0qVzadp9ZGZmVaqkHv0I4F3AJPZ/8Mhnqt1pRLy6l2XzgfnVbtvMsjdv2p3MX3Ve1mEMCJLage8UNR0LfAo4hjL3G/VXJV03t5MMe9wN/K7oZWYtqhlDNfMyHDRrEbE2IqZHxHTgZOBF4FbK3G9UjUq6bsZHxDnV7sDMGqPcw0qamYB9V2zdnQmsj4gngCeK2u8H3l3tRis5o/+ppD+odgdmNjh4xE1ZYyUtLXrN7mXdWcC3S7R/GPh+tQGUPaOX9CjJGPahwIckbSB58IiASP+cMLMG6u/Z+UDvTimcWd34+2bTTjHiyeGVrr6tkhumJA0HZtKji6bE/Ub91lvXja+kmDVZcXdMrUm7Gc+h7b7o6rP5ungb8FBE7CstU+Z+o34rm+jTPiIk/VdE/FnxMkn/BfxZyTeaWU2qTfBZnc07ydfN+yjqtim63+jNPe436rdKLsZOK55J69+cXMtOzaw+BnpXjSUkjQLOBv6yqLnk/UbVbL+3Pvq5wD8AB0vazst1aHaS1pgxs9bWjO4bq116xv7KHm1l7zfqr7KjbiLicxExBvhCRBwaEWPS1ysjourxnGbWXD7rt0qGV/6DpHdK+ldJ/yLpgkYHZWYDn/vuW0clif4rwCXAo8BK4BJJX2loVGZmVjeVXIx9M8ltuAEg6QaSpG9mZgNAJWf0a4GJRfMTgOoebWNmmRhoQzatvio5o38lsFpSd1341wE/k3QHQETMbFRwZjZ4uE+/cSpJ9J9qeBRmZtYwfSb6iPiRpN8DpkTEEkkHA0Mj4reND8/MzGrVZx+9pL8Abga+mjaNB26rZaeSviBpjaQVkm6VdHjaPknSS5KWp69ra9mPmb2sfcLTtE94mhHDTnTf+yBTycXYS4HTgO0AEfEYcGSN++2toP767iL81d7ua2ZmL6sk0XdFxM7uGUlDScoXVy0i7omI3ens/SR/JZiZWQNUcjH2R5K6a96cDfw18N91jOHD7P+8xMmSHib5C+IfI+Inpd6UFu+fDdDW1kahUKhjSPXX2dnpGOvAMdbu5fhupGtX7yOlt6wv1LSv43ZcUPG6hede3lerH8OBppJEPwf4CMlNUn8JfA/4el9vkrQEGFdi0byIuD1dp2dB/S3AxIj4laSTgdskTYuI7T03EhELSIurtbe3R0dHRwX/lOwUCgUcY+0cY+2K41u76cJe1621L3/+qsofJDKraHhlqx/DgaaSUTd7Jd0G3BYRz1W64Yg4q7flpQrqR0QXyVOsiIhlktYDxwNLK92vmZntr2wfvRJXStoGrAHWSnpOUs3j6osK6s8sLqgvqS2td4+kY4EpwIZa92dmNpj1djH24ySjbV6XliY+Avgj4DRJn6hxv18GxpAU1C8eRvkmYIWkR0iGdF4SEc/XuC8zy4jvdm0NvXXdfAA4OyK2dTdExAZJFwH3AF+sdqflCupHxCJgUbXbNbP+8YNJBofezuiHFSf5bmk//bDGhWRmZvXUW6LfWeUyMzPrB0mHS7o5rRiwWtIbipZdLikkja12+7113fxh+qzYA2ICRla7QzMzO8CXgB9ExLslDQdGAUiaQPLQ8Cdr2XjZRB8RQ2rZsJlZpQbzRVtJh5IMRPkgQFqJoLvX5IvA3wO317KPSkogmNkg5MJndTNW0tKi1+wey48FngP+n6SHJX1d0iGSZgJPRcQjtQZQyZ2xZmZWZMhOGPNExSW/tkXEKb0sHwqcBFwWEQ9I+hJwJclZ/ltqCjTlM3ozs2xtBjZHxAPp/M0kiX8y8IikjSSFHx+SVKqsTJ+c6M3MMhQRW4FNktrTpjOBhyLiyIiYFBGTSL4MTkrX7Td33ZiZZe8y4FvpiJsNwIfquXEnejOzjEXEcqBsP356Vl81d92YmeWcE72ZWc450ZuZ5ZwTvdkg5xuj8s+J3sws55zozcxyLpNEL+mzklakT5e6R9LRRcvmSlonaa2kt2YRn5lZnmR1Rv+FiDgxIqYDdwKfApB0AjALmAacA/xH9zNkzcysOpkk+ogornN/CNBdHeh84KaI6IqIx4F1wKnNjs/MLE8yuzNW0nyS59L+BpiRNh8D3F+02ua0zczMqqSIiktt9m/D0hKgVKW1eRFxe9F6c4GREXGFpK8AP4uIb6bLrgO+lz40vOf2ZwOzAdra2k5euHBhI/4ZddPZ2cno0aOzDqNXjrE+Wj3GUvF17Vqx3/yIYSfWbX9bd6zrdfm4ka8+oK2Rx3DGjBnL+igb3KdD2ibE75//iYrWXXbd39W8v1o17Iw+Is6qcNUbgbuAK0jO4CcULRsPlBzkGxELgAUA7e3t0dHRUXWszVAoFHCMtXOMtSsV39pNF+43X8+x9fNXXdPr8lklni7V6sdwoMlq1M2UotmZwJp0+g5glqQRkiYDU4AHmx2f2WDjm6byLas++qvS2st7gSeASwAiYpWkhcAvgN3ApRGxJ6MYzcxyIZNEHxHv6mXZfGB+E8MxM8s13xlrZpZzTvRmZjnnRG9m+/GF2fxxojczy5ikjZIeTet/LS1qvyyt+7VK0tXVbt/PjDUzwGfyLWBGRGzrnpE0g6QszIkR0SXpyGo37DN6M7PW9FfAVRHRBRARz1a7ISd6M7PGGitpadFrdol1ArhH0rKi5ccDp0t6QNKPJL2u2gDcdWNm1k9DdgSHre+qdPVtFdS6OS0ink67ZxZLWkOSn18BvB54HbBQ0rFRRYEyn9GbmWUsIp5Ofz4L3EpSnn0zcEskHiSpJDC2mu070ZuZZUjSIZLGdE8DbwFWArcBZ6TtxwPDgW1lNtMrd92YmWXrVcCtkiDJyTdGxA8kDQeul7QS2AlcXE23TfdGzcwsIxGxAfjDEu07gYvqsQ933ZiZ5ZwTvZlZzjnRm1lm5pV4upTVnxO9mTWUk3n2nOjNzHIuq2fGflbSirRS2z2Sjk7bJ0l6KW1fLunaLOIzM8uTrM7ovxARJ0bEdOBO4FNFy9ZHxPT0dUk24ZmZ5UcmiT4ithfNHkJS0MfMzBogsxumJM0HPgD8BphRtGiypIeB7cA/RsRPyrx/NjAboK2tjUKh0NiAa9TZ2ekY68Ax1i6L+I7bcUHJ9sJzpeNo9WM40DQs0UtaAowrsWheRNweEfOAeZLmAh8FrgC2ABMj4leSTgZukzStx18AAETEAmABQHt7e3R0dDTqn1IXhUIBx1g7x1i7LOKbv+qaku2zyozIafVjONA0LNFHxFkVrnojcBdwRVpgv7vI/jJJ60lqMi/t5f1mZtaLrEbdTCmanQmsSdvbJA1Jp48FpgAbmh+hmVl+ZNVHf5WkdpL6yk8A3aNr3gR8RtJuYA9wSUQ8n1GMZma5kEmij4h3lWlfBCxqcjhmZrnmO2PNzHLOid7MLOec6M3Mcs6J3sws55zozcwyJmmIpIcl3ZnOT5d0f1rccamkU2vZvhO9mVn2PgasLpq/Gvh0WvjxU+l81ZzozcwyJGk88CfA14uaAzg0nT4MeLqWfWRW1MzMzAD4P8DfA2OK2j4O3C3pGpIT8jfWsgMnejOzftKOnQxfs7nS1cdKKq7XtSAtyoik84Bn09peHUXr/BXwiYhYJOm9wHVApfXDDuBEb2aZGETPkt0WEaeUWXYaMFPSucBI4FBJ3wTeTtJvD/Bd9u/W6Tf30ZuZZSQi5kbE+IiYBMwCfhgRF5H0yb85Xe0M4LFa9uMzejOz1vMXwJckDQV2kD5kqVpO9GZmLSAiCkAhnf4f4OR6bdtdN2ZmOedEb2aWc070ZmY5l2mil3S5pJA0tqhtrqR1ktZKemuW8ZmZ5UFmF2MlTQDOBp4sajuBZIjRNOBoYImk4yNiTzZRmpkNfFme0X+R5LbfKGo7H7gpIroi4nFgHVBT1TYzs8EukzN6STOBpyLiEUnFi44B7i+a35y2ldrGbNKxpW1tbRQKhcYEWyednZ2OsQ4cY+2yiO+4HRcc0FZ4rnwMrX4MB5qGJXpJS4BxJRbNA/4BeEupt5VoixJtpLUiFgC0t7dHR0dHdYE2SaFQwDHWzjHWLov45q+65oC2Wb2UQGj1YzjQNCzRR0TJAjyS/gCYDHSfzY8HHkoL628GJhStPp4ay3OamQ12Te+jj4hHI+LIiJiU1nfYDJwUEVuBO4BZkkZImgxMAR5sdoxmZnnSUiUQImKVpIXAL4DdwKUecWNmVpvME316Vl88Px+Yn000Zmb54ztjzcxyzonezCznnOjNzHLOid7MLOec6M2s6QbR82JbghO9mVnOOdGbmeWcE72ZNZW7bQ4kaYikhyXdmc4fIWmxpMfSn6+oZftO9GZm2fsYsLpofg5wb0RMAe5N56vmRG9mliFJ44E/Ab5e1Hw+cEM6fQNwQU37iChZBXhAkfRbYG3WcfRhLLAt6yD64Bjro9VjbPX4oLEx/l5EtNWyAUk/IImxEiOBHUXzC9Iy693buhn4HDAGuDwizpP064g4vGidFyKi6u6bzGvd1MnaiDgl6yB6I2mpY6ydY6xdq8cHrR9jRJxTj+1IOg94NiKWSeqoxzZLyUuiNzMbiE4DZko6l+TM/1BJ3wSekXRURGyRdBTwbC07cR+9mVlGImJuRIxPq/jOAn4YEReRPJvj4nS1i4Hba9lPXhL9gr5XyZxjrA/HWLtWjw8GRoyNdBVwtqTHgLPT+arl4mKsmZmVl5czejMzK8OJ3sws5wZ8opd0jqS1ktZJqunusXqRNEHSfZJWS1ol6WNp+5WSnpK0PH2dm2GMGyU9msaxNG2r623XNcbXXnSclkvaLunjWR9DSddLelbSyqK2ssdN0tz0d3OtpLdmGOMXJK2RtELSrZIOT9snSXqp6Hhem2GMZT/bLI5jrkTEgH0BQ4D1wLHAcOAR4IQWiOso4KR0egzwS+AE4EqSGyJa4dhtBMb2aLsamJNOzwE+n3WcRZ/zVuD3sj6GwJuAk4CVfR239DN/BBgBTE5/V4dkFONbgKHp9OeLYpxUvF7Gx7HkZ5vVcczTa6Cf0Z8KrIuIDRGxE7iJ5NbhTEXEloh4KJ3+LUkNi2Oyjaoidb3tuo7OBNZHxBNZBxIRPwae79Fc7ridD9wUEV0R8TiwjuR3tukxRsQ9EbE7nb0fGN/oOHpT5jiWk8lxzJOBnuiPATYVzW+mxRKqpEnAa4EH0qaPpn8+X59l1wgQwD2Slkmanba9KiK2QPJlBRyZWXT7mwV8u2i+VY5ht3LHrVV/Pz8MfL9ofnJaOfFHkk7PKqhUqc+2VY/jgDHQE71KtLXMeFFJo4FFwMcjYjvwf4HjgOnAFuBfsouO0yLiJOBtwKWS3pRhLGVJGg7MBL6bNrXSMexLy/1+SpoH7Aa+lTZtASZGxGuBvwVulHRoRuGV+2xb7jgONAM90W8GJhTNjweeziiW/UgaRpLkvxURtwBExDMRsSci9gJfI8M/PyPi6fTns8CtaSzPpLdbU4/bruvkbcBDEfEMtNYxLFLuuLXU76eki4HzgPdH2vmddof8Kp1eRtL/fXwW8fXy2bbUcRyIBnqi/zkwRdLk9MxvFsmtw5mSJOA6YHVE/GtR+1FFq70DWNnzvc0g6RBJY7qnSS7UraTOt13Xyfso6rZplWPYQ7njdgcwS9IISZOBKcCDGcSHpHOATwIzI+LFovY2SUPS6WPTGDdkFGO5z7ZljuOAlfXV4FpfwLkko1rWA/OyjieN6Y9J/rRcASxPX+cC/wU8mrbfARyVUXzHkoxieARY1X3cgFeSPOTgsfTnERkfx1HAr4DDitoyPYYkXzpbgF0kZ5of6e24AfPS3821wNsyjHEdST939+/jtem670p/Bx4BHgLenmGMZT/bLI5jnl4ugWBmlnMDvevGzMz64ERvZpZzTvRmZjnnRG9mlnNO9GZmOedEbxVT4n8kva2o7b2SftDkOK6UdHk6/RlJZ1W5nen1qn6ZHpsf9ueuUknnSfp0PfZv1hsneqtYJGNxLwH+VdLI9Gar+cCljdpnmkDL/p5GxKciYkmVm59Ocn9DPZwLPBJJqYtK3UXyYOhRdYrBrCQneuuXiFgJ/DfJXZZXAP8ZEeu7l0t6VVrv/JH09ca0/W8lrUxfHy9a/4D2tEb6akn/QXITzwRJ89Ja5EuA9qL3f0PSu9PpjZI+LekhJbX2p6btp0r6aVq466dKat0PBz4D/Gla+/xP0zuGr5f083Td89P3T5P0YLreCklTShya95PeEZvGv0bS19N/17cknSXpf5XUrD81PZYBFEjKEpg1TtZ3bPk18F7AISR3KD4KjOix7DskRdwgqSN/GHByuu4hwGiSOzFf20v7JGAv8Pp0O93rjQIOJbnL8/J02TeAd6fTG4HL0um/Br6eTh/Ky7XYzwIWpdMfBL5cFPs/Axel04eT3HF9CPDvJPVhIHnuwcEljskTwJh0ehJJ4bA/IDmZWgZcT1Kc63zgtqL3vR/496w/U7/y/Rrar28FMyAififpO0BnRHT1WHwG8IF0vT3AbyT9MXBrRPwOQNItwOkkia9U+x3AExFxf7rN09P1XkzX662e0S3pz2XAO9Ppw4Ab0jPxAIaVee9bSLpSLk/nRwITgZ8B8ySNB26JiMdKvPeISJ490O3xiHg0jXcVcG9EhKRHSb4Iuj0LHN3Lv8esZu66sWrtBfZKmp92aSzvZd1SZWZ7awf4XY/5Smt1dH/x7IF9JzKfBe6LiNcAbydJ4OXieVdETE9fEyNidUTcSFIq+SXgbklnlHjv7h7XEoq/APcWze8tios0lpcq/LeZVcWJ3moSEfO6E2PadC/wVwCShqSjUH4MXCBpVHoB9x3AT3pp7+nHwDskHZxW3Xx7P8M8DHgqnf5gUftvSR712O1u4LK0+iiSXpv+PBbYEBH/RvLXxokl9rGWpFhcfx1Pa1TgtBxzord6+xgwI+2iWAZMi+Sxit8gKS37AEnf+cPl2ntuMF3vOyRVFxdR+sugN1cDn5P0vyTXDbrdB5zQfTGW5Mx/GLBCyUOrP5uu96fAyvSvlqnAf5bYx11ARz/jApiRvtesYVy90qwO0lrq/xkRZ/fjPa8CboyIMxsXmZnP6M3qIpJnxX6tPzdMkVzo/bsGhWS2j8/ozcxyzmf0ZmY550RvZpZzTvRmZjnnRG9mlnNO9GZmOff/AY50D22BvIOcAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.contourf(x_plot, d_pt_arr,trace_sect)\n",
    "plt.colorbar()\n",
    "plt.xlabel('Y-coordinates (m)')\n",
    "plt.ylabel('Depth (m)')\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
